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<£) A technique for the detection and removal of local defects in digital continuous-tone images. 

© The present invention is a method for automatically detecting and correcting a wide range of local digital 
image defects with minimal user intervention. The detection process employs brightness and color thresholds in 
conjunction with magnitude thresholds on residuals of nonlinear spatial filters to separate defects from scene 
content with minimal confusion. The detected defects are then cosmetically corrected by combinations of 
nonlinear smoothing and grey-scale erosion. Several options are outlined for the feature selection, detection, and 
cleaning operations depending on source type and computational constraints. 


/ORIGINAL IMAGE (MUMBANDt) 

V Z 



-~ cleaning block \~180 
(cleaned iuage^, 

FIG. 2 


my onnr.tr>: <ep oe^aasA? i 


Rank Xerox (UK) Business Services 

13. 10/3-09/3.3.41 


EP 0 624 848 A2 


Cross-referenc to Related Applications 

The present application is related to U. S. Patent Application Serial No. 07/934,089 entitled "Process for 
Detecting and Mapping Dirt on The Surface of a Photographic Element" Hied 08#1/92 by Robert Gray et 
5 al.. corres. to EP 0 585 759 A1 . 

Field of the Invention 

The present invention is related to the field of scanning images, primarily from photographic film, to 
w form electrical equivalents of the scanned images and more particularly to the detection and the removal of 
defects due, for example, to dirt contamination or physical damage of the film. 

A portion of the disclosure of this patent document contains material which is subject to copyright 
protection. The copyright owner has no objection to the facsimile reproduction by any one of the patent 
disclosure, as it appears in the Patent and Trademark Office patent files or records, but otherwise reserves 
is all copyright rights whatsoever. 

Background of the Invention 

Digital images created from electronic scanning of continuous-tone photographic film often reveal 

20 visually objectionable defects due to dirt contamination or surface damage of the film. These defects 
usually must be digitally corrected by manual digital retouching (e.g. via "dustbusting" and other cloning 
techniques), which require that the operator visually identify each local defect in the image. 

A patent of interest for its teaching in this art is U.S. Patent No. 4,189,235. entitled "Test Device For 
Dynamically Measuring The Degree of Dirt Accumulation On Bank-Notes" by Guter et al. That patent 

25 describes a method for inspecting opaque web materials for dirt accumulation. A reflected signal from the 
web material is sensed by three adjacent photosensors as the material moves past the sensors. When the 
signal sensed by the center sensor is significantly different than from that sensed by the side sensors the 
sensed area is classified as containing dirt. 

In conventional optical printing systems, if dirt or a scratch appears on the photographic image, the 

30 photographic negative or slide may be cleaned and the print remade or the print itself may be retouched. In 
the case where the images are written to an optical disc, it is not convenient nor desirable to rewrite the 
image, because many images are written onto a disc before the images are viewed, and once an image is 
written it cannot be removed. Therefore, it is desirable to monitor the dirt present in photofinishing 
environments and to assess the effectiveness of any selected film cleaning methods. 

35 In U.S. Patent No. 4,907,156, entitled "Method and System for Enhancement and Detection of Abnormal 
Anatomic Regions in a Digital Image" by K. Doi et al. there is disclosed a technique for taking the 
difference between an "enhanced" and a "degraded" version of the input image to identify and remove 
structured anatomic background. A series of single thresholds are then applied to the difference images 
using the shape/size behavior of image regions to perform a feature extraction that is related to the 

40 detection of lung nodules. The feature extraction is implemented by the use of thresholds against which 
difference image pixels are compared to identify which super threshold "blobs" are lung nodules and which 
are not. 

A significant problem in the detection of dust is to provide an adequate way of detecting "small" 
objects. There are several possible approaches to this problem. For example, one could color-classify all 

45 the color pixels in an image, then spatially connect the pixels of a common color-class into segments, and 
then examine the spatial dimensions of the segment. This approach, while possible, has the drawback that 
the color classification must be quite accurate, or else some marginally-colored pixels will be misclassified. 
The problem is that differences in the direct RGB (or luminance-chrominance) image values do not alone 
give much contrast information between, say, a white small dust artifact and a bright grey extended scene 

so region on which the artifact is superimposed. Thus after color-classification one may be left with very large 
segments that are reasonably neutral and bright, but which now mask visible dust artifacts. 

A better approach prior to classification is to preselect only those points which have large local contrast 
from their surround. This may be directly accomplished by (1) spatially smoothing the color image and (2) 
subtracting the smoothed color image from the original color image. The resulting difference (or "residuals") 

55 image will emphasize all color pixels which have a large magnitude compared to some weighted local 
average magnitude. The key decision here is then the choice of spatial smoother. The most obvious first 
candidate is a low-pass linear filter (e.g. a simple spatially weighted averaging of all grey levels (GLs) in a 
local window); the properties may be easily analyzed mathematically, and (perhaps more important) such 
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filters may be rapidly implemented in software and hardware (the residual image is basically a high-pass- 
filtered version of the original). Unfortunately, linear residuals images have the property of retaining all local 
edges in an image, including very long edges due to scene discontinuities. The result is a very "cluttered" 
image (i.e. the residuals image contains many bright pixels having no connection with dirt or scratches); a 

5 significant amount of additional time must then be spent in scene reasoning to clean up the residuals map 
before defect decisions can be made. In short, linear-filter residuals poorly discriminate between local 
defects and true scene details. 

The identification of anomalies such as dust or scratches leads to the desirability of removing the same 
without the intervention of an operator. The present invention is directed to a methodology for identifying 

io and cosmetically correcting such anomalies with minimal operator intervention. 

Summary of The Invention 

One of the preferred methods of the present invention is: A method for the detection and removal of 
75 defects in digital images comprising the steps of: 

a) preprocessing a digital image represented by pixel values to form a feature image in which defect 
pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of pixel 
values expected, of defects; 

20 c) edge-preserving spatial filtering of each of the feature image pixel values that are within the range of 
expected defects; 

d) forming a residual value for each pixel as a function of the difference between the corresponding 
feature image pixel value and the filtered value of step c); 

e) testing each residual value of step d) to determine if the residual value is within a range of residuals 
25 values expected of defects; 

f) forming a map in which each pixel that is within the range of expected defects according to steps b) 
and e) is marked as a defect; and 

g) correcting the digital image using the map created in step f). 

. The present method, is called an autodustbuster (ADB) algorithm, and it has been successfully tested 
30 on both consumer-type color photographic images and color motion picture frames with good results. 

From the above it can be seen that it is a primary object of the present invention to provide an 
automatic local defect detection and removal technique for digital images. 

The above and other objects of the present invention will become more apparent when taken in 
conjunction with the following description and drawings wherein like characters indicate like parts and which 
35 drawings form a part of the present invention. 

Brief Description of The Drawings 

Figure 1 illustrates a system on which the present invention may be implemented. 
40 Figure 2 illustrates the windowing of an image raster to determine the grey level of a central pixel GLc 
at a position X. 

Figure 3 illustrates in chart form a histogram of GL values for an associated window. 
Figure 4 illustrates the averaging of the samples "nearest" in grey level GL to the central pixel GLc for 
the histogram of Figure 3. 

45 Figure 5 illustrates the replacement of the GLc value with the grey level average GL av at the position X. 
Figure 6 illustrates the process flow for a multi-band original image. 
Figure 7 illustrates the process flow for a resolution pyramid of an image. 
Figure 8 illustrates the process flow as applied to an R, G, B image. 
Figure 9 illustrates a flow for ADB processing within a resolution pyramid. 

so 

Detailed Description of the Invention 

Referring to Figure 1, the present invention is implemented using a digitizing scanner 10, a control and 
logic device 24, such as a computer with a display device 26, and an optical disc writer 36. As shown, a 
55 strip of film 15, containing one or more frames 18 of developed film, is placed into a scan gate 12 where 
the strip of film 15 is scanned under direction of control and logic circuitry 30. As each frame 18 is scanned 
the resultant scanned image data, represented by block 22, is digitized and transmitted to a memory 28 for 
storage. The computer 24 process the stored data, in a manner to be described, to provide output image 
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data 32 which may be written to an optical disc 34 by the optical disc writer 36 to provide a report as to the 
characteristics of the anomalies. The scanning device 10 is capable of quantizing pixel values into multiple 
brightness levels in separate red, green, and blue channels. A minimum number of brightness levels would 
be approximately 64 for an adequate quality image with a typical number being 256. Although the image 

5 data 22 has been described as originating from the operation of scanning film, it is well understood that 
other well-known techniques may be used to provide the image data. 

The autodustbuster (ADB) algorithm and several processing options that allow customization for 
different source types and computational constraints will now be described. 

The required ADB operations may be divided into (1) feature selection, (2) defect detection, and (3) 

io defect correction. The operations must ensure that the maximum number of defects is detected and 
optimally corrected while minimally affecting scene elements which are not defects. Thus, the more 
confidence which the user has in the accuracy of the detection process, the more extreme can be the 
removal process; conversely, as one's confidence in avoiding "false alarms" from scene elements 
diminishes, milder removal methods may be used to avoid damaging real scene elements. 

15 Visible dirt artifacts in digital positives, created from scanned color negatives, possesses the following 
properties: they are "white" (low color saturation), "bright" (very high apparent exposure values), "small" 
(small spatial extent) and "contrast" (large differences from their immediate surround). {If dirty color 
positives were instead scanned, then the dirt artifacts would be "dark" rather than "bright"). In addition, 
many such fragments exhibit little brightness variation within individual defects. 

20 It is obvious that the above characteristics are rather loosely defined, and that many real scene 
characteristics may possess one or more of these attributes. Thus the goals of an ADB algorithm are to 
maximize the percentage of defects detected ("hits") while minimizing the number of real scene elements 
mistakenly considered to be defects ("false alarms"). These goals tend to be competitive, that is, setting 
thresholds on features to increase the percentage of "hits" tends to increase the number of "false alarms" 

25 and vice versa. A detection mechanism must be used in which defects appear as dissimilar as possible 
from true scene elements. 

In most imagery, spatially small objects which differ greatly in brightness or color from the immediate 
surround are relatively rare. The ADB algorithm exploits this property by the innovative use of a particular 
form of nonlinear spatial filter. The properties of this filter are described in the next section. 

30 Figure 2, illustrates in flow chart form the generic autodustbuster process flow. The original image of 
block 100 is a multi-band digital image represented by one or more matrices of pixel values or equivalents 
grey levels (GL). In preprocessing block 110 some effective and convenient remapping of the original image 
may be performed. An example of preprocessing is the digital resampling of the image to reduce its size 
for computational expediency. An alternate preprocessing is the selection or recombination of the original 

35 image bands to produce a feature image with fewer bands. The detection suitability of the resulting feature 
image, in any case, should be at least comparable to that of the original image. Oval 120 represents the 
feature image that is outputted from block 110. The feature image is subjected to a band brightness test in 
block 130. The test applies low and high threshold evaluations to every pixel of each band in the feature 
image. Pixel values that are above the high threshold or below the low threshold in each band are 

40 considered to have passed the test and are forwarded to an edge-preserving-smoothing (EPS) filter 140. 
This filter non-linearly smoothes the selected pixel in all bands in a manner that retains major scene edges 
as subsequently described. The output of the filter at this pixel position is written to a corresponding 
position in the smoothed image 142. For those pixel values that are below the high threshold value, but 
above the low threshold value no subsequent EPS filtering is applied, i.e. the test is considered failed. Each 

45 band of the feature image pixel is copied to the resultant output of the smoothed image 142 and does not 
undergo EPS filtering. At difference node 145 the pixel values of each band of the smoothed image are 
subtracted from the pixel values at the corresponding position in the feature image 120 resulting in the 
residuals image 148. 

All of the pixel values from the residuals image are subjected to the residuals test of block 150. The 
so residuals test consists of a minimum and a maximum threshold that is applied to each band of the residuals 
image 148. When bright defects are to be detected, for a pixel to pass all band values must exceed the 
minimum threshold and at least one band value must exceed the maximum threshold. When a dark defect 
is to be detected a pixel will pass if all band values are less than the maximum threshold and at least one 
band value is less than the minimum threshold. 
55 Pixels that pass are forwarded to a classmap 160 and marked with a distinctive value. The classmap 
160 is a single band image representation which is used to mark the location of detected defects. Pixels 
that fail are marked with a different value in the classmap 160. In a classmap conditioning block 170 the 
values within the classmap 160 may be further modified to accommodate known spatial characteristics of 
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the defects. For example, a one-pixel morphological dilation of the defect-marked (passed) pixels in the 
classmap is typically performed. The conditioned classmap 175 is directed to the cleaning block 180 long 
with the original image pixel values. Pixels in the original image that are not marked as passed in the 
classmap will be copied unchanged to the output cleaned image 190. Pixels which are marked as passed, 
i.e., defect, are corrected using neighboring non-defect pixel values. The corrected values are then written 
to the cleaned image 190. 

A more detailed discussion of the major blocks of Figure 2 will now be undertaken. In ADB processing 
the feature image 120 used to detect the defect locations may be different from the original image 100 
which will be corrected. This is generally advantageous when computational time is at a premium and EPS 
filtering 140 of the full (generally three-color) image is too burdensome. In this case it is often possible to 
apply a preprocessing module 110 to create a one- or two-band "feature" image 120 which instead will be 
EPS-filtered and classified for defects. The effectiveness of this strategy will depend on the success with 
which a reduced-dimensionality feature image will retain the size, brightness, contrast, and color defect 
discrimination of the original image 100. This will of course be application-specific. For example, in digitized 
imagery originating from superimposed color cells (as in animated motion pictures) defects may be of any 
color and saturation; thus little is lost in detection performance using a reduced-dimension space jwhich 
sacrifices color and saturation information. In this case good success has been achieved in using a "band- 
maximum" one-band feature image 120 band where every feature pixel value is the maximum value of the 
input red. green, and blue bands. By contrast, in Photo CD applications the principle source of defects is 
dirt on the negative during scanning, which appears as "white" specks in the digital positive; here saturation 
information should be retained as an important identifier, and two- or three-band EPS-filtered feature images 
are the norm. Figure 7 illustrates the retention of all three color original color bands for classification. 

Rescaling of the feature image grey-level range during preprocessing 110 is often advantageous. In 
current implementations this consists of linearly scaling to a 0-255 pixel value range even for 10- or 12-bit 
original data. This results in a slight speed improvement during EPS-filtering 140 when the histogram is 
searched at every window location (see Figures 4 and 5 along with associated discussions). In Photo CD 
applications (e.g. Figure 7 this rescaling may be determined by finding .5% and 99.5% cumulative 
histogram points of the color-corrected bands of the entire input image. The minimum .5%-point and 
maximum 99.5%-point from among all bands are mapped to 0 and 255 pixel values, respectively in the 
feature image, thus ensuring that no color shift occurs. If the input data is not color balanced it may be 
necessary to perform a color balancing step on the original data. 

EPS Filters and k-Nearest Neighbor (KNN) 

To detect the critical size of the desired defects a smoothing filter is desired which retains extended 
image edges but smoothes small local contrasty image regions (blobs). Such filters exist as various 
nonlinear noise-reduction methods under the category of "edge-preserving smoothers" (EPS). The most 
widely known EPS filter is the median filter, which is one member of the class of rank-order filters. In rank- 
order filters, the pixel magnitudes are sorted and an average of some contiguous number of the sorted 
magnitudes is computed. For example, in a median filter the "middle" element of the sorted magnitudes is 
retained (an average over one sample). The sample mean is, at the other extreme, also a rank-order filter 
(although a poor EPS); it provides the average of all sorted magnitude values. 

A poor property of standard rank-order filters for image operations is that they are insensitive to the 
relative spatial location of the pixels within the current window; e.g. the magnitude chosen as the "median" 
may exist at the center of the window or at an extreme edge of the window. A number of specialized 
nonlinear window operators have-been devised to account for local spatial structure. One example is the 
sigma filter (J.-S. Lee, CVGIP.24,255-269. 1983). In this filter a histogram is created of the grey levels in a 
sliding window; based upon an estimate of the image noise, a range of histogram values about the grey 
level of the pixel at the center of the window is computed. One characteristic of the sigma filter, however, is 
that it is specifically designed for noise removal and its averaging range is solely dependent on the noise 
characteristics; the filter averages over a different number of pixels at each window location. For "blob- 
detection, however, what is desired is a filter whose parameters are sensitive to local shape characteristics. 

A relative of the sigma filter which gives this control is the k-nearest-neighbor (kNN) filter (Davis and 
Rosenfeld, IEEE Trans, SMC, 7, 107-109 ,1978). This filter is implemented in the EPS filter module 140 of 
Figure 2. This filter has the following operation per color band: 

Referring to Figure 3. slide a rectangular window of odd dimensions n x m over an image in raster 
order. At each location determine the grey level of the center pixel (GU) or equivalents pixel value, and 
also create a histogram of the grey levels in the local window as per Figure 4. -Beginning at the histogram 
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bin containing the center grey level GL Cl sum contents of the histogram bins on both sides of the starting 
bin in the following fashion: Test if the, contents in the bin for Gl_c are greater than or equal to k pixels. If 
this is true, then write out as the averaged value GL av just the existing value GU. If the contents for bin Gl_c 
are less than k, then add the contents of the bin immediately below the bin and above Gl_c to the sum and, 
s again test whether the sum is greater than or equal to k. If the sum equals k, then calculate the average 
grey level according to the formula 

GL av = SUM {GLi x hist(ij) 

w 

SUM {hist (i) } 


75 where the index i is over the summed histogram bins. If the sum is greater than k, then the extreme bins in 
the sum are proportionally reduced in weight in the calculation for GL av so that the total bin weight in the 
sum approximately equals k (as shown in the attached source code). If the sum is less than k, then the next 
bin below GLc and above GU are summmed and tested against k. This process is repeated until the sum of 
histogram bins about GL C are greater than or equal to k. Continue until k pixels <k < n x m) have been 

20 summed. Find the average grey level of the summed pixels per Figure 5 and write this value to an output 
image at the location of the window center per Figure 6. 

The k-nearest neighbor filter has some interesting filtering properties which are not immediately obvious 
from the above description. First, any constant-GL image structure (blob) which has an area of k or more 
pixels in the local window will not be significantly smoothed: any structure of less than k pixels will be 

25 increasingly smoothed as k decreases. Second, at "ramp" edges between otherwise constant grey level 
regions the slope of the ramp will be increased (sharpened) by the filter. 

The first property is the one of significance for "blob" detection. For example, if the window size is 5 x 
5 pixels, and if k is, say, 15, then extended straight edges having a depth (constant pixel value in both 
directions perpendicular to the edge) greater than three pixels will be effectively unchanged, and thus 

30 disappear in a residuals image. "Blobs" of area less than 15 pixels that fit entirely in the window will be 
smoothed and thus will appear strongly in the residuals image. This is the desired residuals property for a 
local defect detector. 

A major concern for a k-nearest neighbor filtering is computational speed. This is not a linear filter and it 
is not separable, nor does it have other usable symmetry properties. The major operations per window 

35 location are histogram creation and then accumulation with testing of the total accumulated samples. The 
current k-nearest neighbor implementation does significantly reduce the cost of computing the histogram by 
merely updating the window histogram in a horizontal strip of image rather than performing complete 
calculations per window. This results in a fractional saving of approximately (n-2)/n in histogram-computing 
operations (for a n pixel x m line window), i.e. the savings increase as window length increases. (This fast 

40 histogram updating is equivalent to that used in the "fast-median" filter-Huang and Tang, IEEE Trans, 
ASSP, 27, 13-18, 1979). 

Residuals Test 

45 As mentioned above, a residuals image 148 is created as shown in Figure 2. The anomaly classifier 
steps contained in modules 150 and 170 accept as input the residuals image 148 and create as output the 
final classmap 175 which marks the pixels to be cleaned. 

The residuals test in module 150 creates an initial classmap by performing tests on each pixel of the 
residuals image 148. In the classmap conditioning module 170 a spatial segmentation (i.e. connected- 

50 component labeling) of the non-background pixels in the classmap may be optionally performed with 
additional tests applied to the same segment. 

In the simplest (single-state) version and for detecting bright defects, each color residual pixel is 
evaluated against two residual-GL thresholds, GL rmin and GL rmax . In this implementation, all color residual 
pixels whose minimum band residual GL exceeds GL rmin andwhose maximum band residual GL exceeds 

55 GWmax are uniquely marked for subsequent cleaning; pixels which do not satisfy both of these constraints 
are marked to be left unaltered in the output image. (For single-band feature images Ql rmm equals GUmax)- 
The two-dimensional map of category markings per pixel is termed the "class map". These thresholds 
impose a minimum required contrast between a defect pixel and its immediate background. The (smaller) 
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GL rmin threshold is generally set slightly above the "zero point" in the residuals image (e.g. to 133 for 
rescaled features with residuals biased by a value of 128 as described above) and acknowledges that valid 
bright defect pixels , e.g. due to dirt on a scanned negative, should not have a negative contrast (residual 
values less than the zero point) in any band. The (larger) GL^ threshold states a requirement of 

5 significant local contrast in at least one band for the defect to be visible. (Recalling that these thresholds are 
on the residuals GLs, it is obvious that not all bands should be required to pass a large GL^ax value; for 
example, a white dirt artifact against a bright red background would have large residual values in the green 
and blue bands, but only a small residual value in the red band.) 

For two-state anomaly classification, two sets of residuals thresholds are used, i.e. GL rmin i, GUnaxi, 

io GL rmin2 , and GL rmax2 . In general, the minimum-residuals thresholds are set equally, and slightly above the 
zero point as mentioned above; the maximum-residuals threshold is set somewhat higher for labeling high- 
confidence pixels than is the maximum-residuals threshold for the low-confidence pixels. This labeling 
results in a class map where the pixels possess one of three labels (normal, low-confidence dirt, and high- 
confidence dirt). In two-state classification, this class map undergoes image segmentation in module 170, 

75 i.e. all residuals pixels spatially contiguous in the class map which have the same class are, assigned a 
unique code value. (These code values can be thought of as uniquely labeling all the "blobs" of suspected 
defects in the image.) A fairly efficient algorithm for segmenting the class map is employed which involves 
finding the lowest prototype code value in a spatial "tree" of prototype segments of the same class, and 
then renaming all the segments in the tree with this lowest code value. 

20 Following segmentation, all segments which contain low-confidence defect pixels and which are 
spatially touching segments containing high-confidence pixels are upgraded to high-confidence (i.e. the 
class map values of the pixels are reassigned to the "high confidence of an defect" class). The purpose of 
this step is to make use of the adjacent presence of very obvious defects to decide that a slightly less 
contrasty region is part of the same defect. The segments which remain at low confidence are either 

25 relabeled as "clean" (i.e. all pixels of the segment are marked for no cleaning in the classmap) or else 
name their pixels marked for subsequent cleaning state in module 180. The segments, which are marked as 
high confidence of being defects, have their pixels marked in the classmap for subsequent erosion cleaning 
in module 180. 

30 Cleaning 

Two general cleaning methods are presently available for use within module 180 in the ADB algorithms: 
EPS substitution and greyscale erosion. In simple EPS substitution, the original input multi-band grey levels 
of image 100 are simply replaced by the band EPS-filtered grey levels of smoothed image 142 wherever a 

35 pixel is marked as a defect in the classmap 175 (as shown in Figure 2). This method has the advantage of 
computational simplicity in that EPS values have already been computed during defect detection; however it 
requires that the EPS image be stored for use during cleaning. In addition, EPS substitution is relatively 
forgiving of "false alarms", i.e. pixels classified as dirt artifacts but which are really part of the scene. 
Conversely, EPS-substitution cleaning may not completely remove all of an extended defect, to the degree 

40 that the defect itself influenced the local EPS values. 

In the greyscale erosion cleaning option, an iterative process occurs over detected defect regions, 
whereby GLs of defect pixels are replaced by a weighted average of original GLs of "normal" (i.e. non- 
defect or corrected defect) pixels within a small distance of the defect pixel. As the iterations progress, the 
effect of the "normal" GLs propagates toward the center of the original extended defect. The iterations end 

45 when all marked defect pixels have been replaced. Two precautions are in order. First, the "normar-labeled 
pixels that are in immediate contact with the defect pixels at the beginning of the iterations are also 
relabeled to class "defect" prior to the onset of cleaning (via a morphological dilation operation which uses 
a 3x3-pixel probe function). This is because such pixels, while failing the detection thresholds, often have 
some contribution from the defect in their GLs which will be visible when propagated over multiple pixels 

so during erosion. Cleaning their values as well as those of "true" defect pixels can significantly improve the 
quality of the cleaning. The second precaution is that the spatial window over which the erosion weights are 
taken should be kept reasonably small, e.g. 5 x 5 samples. The reason for this is to minimize the influence 
of distant pixels which may not be part of the same scene structure as the region obscured by the central 
defect pixel. 

55 Current practice for the erosion cleaner is to use a 5x5-pixel sliding window with a minimum of six 
clean-class contained pixels in order for a central defect-class pixel to be filtered; if the central pixel is 
defect-class, but fewer than six clean samples are contained, then the defect classification is maintained to 
the next iteration. If the central pixel is of defect class and six or more clean-class samples are contained in 
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the window, then the GLs of the central pixel are replaced by a weighted average of the band GLs of the 
contained clean pixels, where the weights are linearly proportional to the inverse of the Euclidean spatial 
distance of the clean pixel from the central pixel. The cleaning is iterative (multiple passes through the class 
map) but is not recursive, i.e. the cleaning effect at one window location does not affect the cleaning at 

s subsequent locations during the same iteration, in order to avoid directional biasing. 

Greyscale erosion has the advantage over EPS-substitution in being able to completely eliminate all 
trace of the artifact from the image; its disadvantages are increased complexity, a variable processing rate 
(the number of iterations is not known a priori), and the possibility for severe scene damage in the case of 
"false alarms" (i.e., real scene elements may be removed from the image). Despite these concerns, the 

70 excellent cleaning performance of erosion is currently the preferred method in both Photo CD (Figure 7) 
and Cineon (Figure 8) applications. 

Alternative Implementations: 

75 The ADB algorithm has several "flavors". These include (a) operation at a single spatial resolution, (b) 
operation within a resolution pyramid, (c) single-defect-state detection, (d) multiple-defect- state detection, 
and cleaning via (e) EPS substitution or (f) grey scale erosion. In the single-state versions, all pixels of the 
image are classified as either "non-defect" or "defect", i.e. only one defect state is treated. In multi-state 
versions, the "defect" classification is further divided into subclassifications depending on the confidence 
20 that a pixel is indeed a defect; thus a three-state classification would consist of rt non-defect M . "high- 
confidence defect", and "transition-confidence defect". In single-state versions the cleaner used is either 
EPS-substitution or grey-scale erosion, but not both in one implementation. The multi-state version may 
employ a combination of both EPS-substitution for cleaning the lower-confidence defect regions and grey- 
scale erosion for cleaning the higher-confidence defect regions. 

The flavor combining options (a), (c). and (f) is a simple variant to implement, though it is not 
necessarily the fastest or best-performing version. In this case no spatial resolution pyramid is formed, the 
residuals image is not segmented, and the cleaner is chosen to be gray scale erosion for maximum 
correction. Figure 7 illustrates such a process flow and is typical of a Photo CD application. In this case the 
preprocessing 1 10 consists of rescaling the pixel values of the original image for maximum viewing contrast 
by linearly remapping the lowest .5-% band gray level to a value of zero and the highest 99.5% level to 
255. as previously described in the description for Figure 2. The brightness test 130 for dirt on scanned 
negatives consists of minimum threshold pixel values which must be equaled or exceeded in each band of 
the feature image 120 in order for the corresponding pixel to pass the test as a possible defect. The 
residuals classification 150 is a single state decision per pixel as either a non-defect or defect class, based 
35 upon a requirement that the smallest band residual equal or exceed the zero-background residuals level 
and that the largest band residual equal or exceed some user-determined threshold level which is above the 
zero background level. The classmap conditioning 170 in this case consists of a morphological dilation of 
the defect pixels in the classmap using a 3 x 3 or 5 x 5 kernel. The cleaning 180 consists of greyscale 
erosion as previously described under the section called "cleaning." 

Figure 8 illustrates an alternate processing flow which retains the characteristics, (a), (c) and (f), is 
appropriate for removal of bright defects in Cineon applications. In this case it is desirable to minimize the 
processing time required by the large image format the processing is therefor directed at only detecting the 
largest defects and at minimizing the number of bands in the feature image. The preprocessing 110 in this 
instance consists of creating a reduced resolution image rescaling it to a zero to 255 greylevel range such 
that the color balance is corrected, and then extracting a single band feature image in which the output 
pixel value is the larger of the Red. Green and Blue input values. The residuals test 150 applies a single 
threshold which each residuals pixel must equal or exceed to be marked as defect. Classmap conditioning 
170 consists of morphological dilation as described above. The cleaning module 180 consists initially of 
performing greyscale erosion of the reduced resolution unrescaled input image as previously descnbed. 
Next the cleaned reduced resolution multi-band image is bilinearly interpolated back to the original size and 
the conditioned classmap is replicated to the original image size. Finally those pixels in the original image 
which are marked "defect" in the replicated classmap are replaced by the pixel values in the interpolated 
cleaned reduced resolution image. 

It is possible to use ADB within a multi-resolution pyramid structure to provide some speed gams from 
use of multiple smaller EPS kernels. This implementation would also allow more effective cleaning from the 
use of EPS residuals thresholds which are tuned to the specific pyramid level. 

Figure 9 illustrates a flow for ADB processing within a resolution pyramid. In most respects. ADB 
cleaning within a spatial level is similar to that in a single-level ADB version, i.e. detection, classification. 
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and cleaning occur, and choices of one- or two-state and of EPS-substitution and/or greyscale erosion are 
the same. One key difference lies in the way the spatial pyramid residuals images are processed. 

In a spatial pyramid, an image is represented generally as a low-resolution image plus a number of 
residuals images of increasing resolution. To reconstruct the image to a particular resolution, one interpo- 

5 lates the low-resolution image by a prescribed method to the next resolution level and then adds the 
residuals values of that level to each pixel GL. This process continues until the desired resolution level has 
been reconstructed. The Photo CD file format is an example of a spatial pyramid structure. Note that these 
residual images are simply the pixel band difference between the original image at that resolution and the 
interpolated low-resolution version; they should not be confused with the EPS residual images discussed 

10 above. 

In a pyramid structure, all higher-resolution reconstructed pixels depend on the lower-resolution GLs at 
that spatial location. Thus when some subset of defect pixels is cleaned at a given resolution level, the 
residuals values at all higher resolutions at these locations must also be adjusted; failure to correct the 
residuals* values could result in "ghost" images of the original defects upon pyramid reconstruction. 

75 Although various residuals processing schemes are possible, the method currently used for ADB is the 
following relatively simple one: the base (lowest level) version of the original image is cleaned as in Figure 
7. The classmap 1 75 prior to base cleaning is then replicated in module 200 to the size of the next highest 
resolution level. The defect-class pixels in the replicated classmap are then dilated in module 210 by a 
spatial extent equal to the region of support of the pyramid smoothing function for that level, i.e., the half- 

20 width of the kernal used to create the next-lower-resolution level. The resulting dilated classmap 220 is then 
compared with the pyramid-residuals image 230; the pyramid residuals image is assigned a zero pixel 
value in every band wherever the corresponding pixel in the class map 220 is labeled as defect. The 
resulting modified pyramid-residuals image 250 is then added per-pixel in module 270 to a version of the 
base level cleaned image 190 which has been spatially interpolated to the present pyramid resolution level 

25 in module 260. * 

This results in a reassembled (base + 1)-level image 100' in which the results of the base-level ADB 
learning have been incorporated/This image is then ADB-processed starting with the preprocessing step 
110\ and the process of Figure 9 thus repeats until all levels have been reassembled and cleaned. Note 
that ADB cleaning may be omitted from higher-resolution levels at the expense of small defects not being 

30 corrected, but that all higher-level, pyramid-residuals must be zeroed as described in order to ensure that 
no defect "ghosts'* will appear in the final image. 

While there has been shown what are considered to be the preferred embodiments of the invention, it 
will be manifest that many changes and modifications may be made therein without departing from the 
essential spirit of the invention. It is intended, therefore, in the annexed claims, to cover all such changes 

35 and modifications as may fall within the true scope of the invention. 
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Scanner device 

40 

12 

Scan gate 


15 

Film strip 


16 

Film cleaner 


18 

Frame 


22 

Image data 

45 

24 

Control and logic device (computer) 


26 

Display device 


28 

Memory 


30 

Logic circuitry 


32 

Output image data 

50 

34 

Optical disc 


36 

Disc writer 


100 

Original image 


no 

Preprocessing 


120 

Feature image 

55 

130 

Band brightness test 


140 

Edge-preserving-smoothing (EPS) filter 


142 

Smoothed image 


145 

Difference or subtracting node 
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148 

Residual image 

150 

Residuals test 

160 

Classmap 

170 

Conditioning block 

175 

Conditioned classmap 

180 

Cleaning block 

190 

Cleaned image 

200 

Classmap replication block 

210 

Classmap dilation block 

220 

Replicated dilated classmap 

230 

Higher-resolution spatial resolution-residuals image 

240 

Spatial residuals zeroing block 

250 

Modified higher-resolution spatial resolution-residuals image 

260 

Cleaned-image interpolation block 

270 

Pixel value summing node 

100' 

Rebuilt modified higher-resolution image 

110* 

Higher-resolution preprocessing 


The Invention may be summarized as follows: 

1. A method for the detection and removal of defects in digital images comprising the steps of: 

a) preprocessing a digital image represented by pixel values to form a feature image in which defect 
pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of 
pixel values expected of defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values that are within the range 
of expected defects; 

d) forming a residual value for each pixel as a function of the difference between the corresponding 
feature image pixel value and the filtered value of step c); 

e) testing each residual value of step d) to determine if the residual value is within a range of 
residuals values expected of defects; 

f) forming a map in which each pixel that is within the range of expected defects according to steps 
b) and e) is marked as a defect; and 

g) correcting the digital image using the map created in step f). 

2. A method for the detection and removal of defects in digital images comprising the steps of: 

a) preprocessing a digital image represented by pixel values to form a feature image in which defect 
pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of 
pixel values expected of defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values that are within the range 
of expected defects; 

d) forming a residuals image which consists of the difference between the feature image pixel values 
and the filtered values of step c); 

e) testing each pixel value of the residuals image of step d) to determine if the residual value is within 
a range of residuals values expected of defects; 

f) forming a map in which each pixel that is within the range of expected defects according to steps 
b) and e) is marked as a defect; and 

g) correcting the digital image using the map created in step f). 

3. A method for the detection and removal of defects in digital images comprising the steps of: 

a) preprocessing a digital image represented by pixel values to form a feature image in which defect 
pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of 
pixel values expected of defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values; 

d) forming a residual value for each pixel as a function of the difference between the corresponding 
feature image pixel value and the filtered value of step c); 

e) testing each residual value of step d) to determine if the residual value is within a range of 
residuals values expected of defects; 
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f) forming a map in which each pixel that is within the range of expected defects according to steps 
b) and e) is marked as a defect; and 

g) correcting the digital image using the map created in step f). 

4. The method wherein the map formed in step f) is conditioned by morphological dilation o( the pixels 
marked as defect. 

5. The method wherein the map formed in step f) is conditioned by morphological erosion of the pixels 
not marked as defect. 

6. The method wherein the map of step f) is formed with each pixel marked according to the magnitude 
of its values according to steps b) and e). 

7. The method wherein the map formed in step f) is formed with each pixel marked according to the 
magnitude of its values according to steps b) and e)and is conditioned by morphological dilation of the 
pixels marked as defect. 

8. The method wherein the map formed in step f) is formed with each pixel marked according to the 
magnitude of its values according to steps b) and e) and is conditioned by morphological erosion of the 
pixels not marked as defect. 

9. The method wherein the filtering of step c) is performed with a k-nearest neighbor filter. 

10. A method for the detection and removal of defects in digital images represented by multiple spatial 
resolutions ordered from lower to higher, comprising the steps of: 

a) preprocessing a lower spatial resolution representation of the digital image by forming a feature 
image in which defect pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of 
pixel values expected of defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values that are within the range 
of expected defects; 

d) forming a residual value for each pixel as a function of the difference between the corresponding 
feature image pixel value and the filtered value of step c); 

e) testing each residual value of step d) to determine if the residual value is within a range of 
residuals values expected of defects; 

f) forming a map in which each pixel that is within the range of expected defects according to steps 
b) and e) is marked as a defect; 

g) correcting the lower resolution representation of the digital image using the map created in step f); 
and 

h) correcting the higher resolution representations of the digital image using the map created in step 
f). 

11. The method wherein one or more of the higher resolution representations previously corrected using 
the lower resolution defect map of step f) are again processed using steps a) through h) to detect and 
correct for defects too small to be detected at lower spatial resolution representations. 

12. The method wherein the filtering of step c) is performed with a k-nearest neighbor filter. 

13. The method wherein the filtering of step c) is performed with a k-nearest neighbor filter. 

14. An apparatus for the detection and removal of defects in digital images comprising: 

means for preprocessing a digital image represented by pixel values to form a feature image in 
which defect pixels have large local contrast from their neighboring non-defect pixels; 

first means for testing the value of each feature image pixel value to determine if the value is within 
a range of pixel values expected of defects; 

filter means for edge-preserving spatial filtering of each of the feature image pixel values; 

means for forming a residual value for each pixel as a function of the difference between the 
corresponding feature image pixel value and the 'filtered value from said filter means; 

second means for testing each residual value from said forming means to determine if the residual 
value is within a range of residuals values expected of defects; 

means for forming a map in which each pixel that is within the range of expected defects according 
to said first and said second testing means is marked as a defect; and 

means for correcting the digital image as a function of the formed map. 

15. An apparatus for the detection and removal of defects in digital images represented by multiple 
spatial resolutions ordered from lower to higher, comprising: 

means for preprocessing a lower spatial resolution representation of the digital image by forming a 
feature image in which defect pixels have large local contrast from their neighboring non-defect pixels; 

first testing means for testing the value of each feature image pixel value to determine if the value is 
within a range of pixel values expected of defects; 
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filter means for edge-preserving spatial filtering of each of the feature image pixel values that are 
within the range of expected defects; 

residual means for forming a residual value for each pixel as a function of the difference between the 
corresponding feature image pixel value and the filtered value of said filter means; 

second testing means for testing each residual value from said residual means to determine if the 
residual value is within a range of residuals values expected of defects; 

mapping means for forming a map in which each pixel that is within the range of expected defects 
from said first and said second testing means is marked as a defect; 

means for correcting the lower resolution representation of the digital image using the map created 
by said mapping means; and 

means for correcting the higher resolution representations of the digital image using the map created 
by said mapping means. 

16. A method for the detection and removal of local defects in digital images comprising the steps of: 

a) generating a feature image from an original digital image; 

b) creating an EPS-residuals image using the feature image; 

c) testing and mapping image pixels based upon the values of the EPS residuals image and the pixel 
values of the feature image; 

d) modifying the mapped image as a function of mapped values of neighboring pixels; 

e) cleaning of defect pixels by an EPS-substitution process or by a greyscale erosion process; and 

f) cleaning of defect pixels in different spatial resolution representations of the original image if 
different spatial resolution representations are present. 
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- Busteps 042193: (4-21-93) Applies a selective-averaging 

(Tc-NN smoothing) filter to a 3-band image; also <„i„ a i 
i calculates the difference (residual image) between the original 

: & the smoothed images - 

Search window size must be symmetric, odd of dimensions 

3x3 thru 51x51. n „ oori . , 

Filter currently treats data as in range 0-4 999 *oniy - 
Input image array is stored as Int*2. 

intended for operation on a Sun workstation (e.g. Sparcstation 10) 
running Sun OS. 


Calls: 
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Read Hdr isl 

— Timage-f ormat-specif ic calls) 

Ihister 

( image -format -specif ic caxls) 

Ihismnmx 
Winbust 042193 

"(image-format-specific calls) 

Hist Init 

Hist^UP date 
Iselavh 

Add Image 

~ (image-format-specific calls) 

program BustepS_042193 
integer*2 npix,nlin 

integer*2 ibuf f (4096, 3) , buf f n (4 096, 51 , 3) 
integer nav, npp, diagf lg 

Integer npixm,nlm, np, nl,band, nbands, nbx 
integer pixf ormi, pixf ormo, pixf ormro, dumml, dumm2 
integer gmin (3) , gmax (3) 
integer hist (5000, 3) 
integer*4 britelo, britehi 

real*4 f actl, f act2, bias 

real*4 deltime, Dtirae, timearray (2) 

character*50 infil, outf il, resf il 

data npm, nlm/51,51/ !Max EPS window size! # « i 

!data pixformro/3/ !residuals image written as unsigned byte! 

data factl, fact2/l. 0,-1-0/ 'Used m Aad_Image! 

data npixm, nbx/4096, 3/ 

data diagf lg/0/ !For winbust3! 

common /pa r/ dumml, dumm2, nav 


1 
2 
3 
4 

5 


format ( 


Interactive parameter input 

write (6,5) 

-Busteps 042193":'/, ^ ^ w . , . . 

Tnput image may be 1,2, or 3 bands;/. 
Input values should be in range 0-4999; /, 
Max. image width: 4096 pixels per line; /, 
Max. window height: 51 pixels'/, 
Left- & right-edge virtual mirroring /) 

write(6,10) ^ crt „ , 

format ( ' Enter in f ilenaroe [ch50] : , s ) ^. . . 

read(5,' (aS0)')infil !Must be ISL IPP imagefile! 

write (6,30) 
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75 


EPS filename [ch50] |:' * s > 
35 read(5,' (a50)') restiJ- >=3 ,«<51) : ' . 5) 

saiva-* 'nafsusas 


35 


40 


Call Exitd) 

Endif , no 

nl=np 

io npp=np*ni 


50 

formate 
51 read(S,«> nav 


write (6 5 50) of avera ging points in window') 

formatt' Enter no. 01 , n t 3 *nl)/2 ]:'.5> 

write (6, 51) x f typically about (np nl)/ 

formate 1 r 


stop 
Endif 


SSJlSl E Srit.lo.brit*hi 

f^' 6 ( :^Ur cut residuals bias constant S ) 
forrnatl Enter .^-..oixfon 


20 

60 

stop |5 

Endif ^ a «t*<; ■ 

initialize arrays * consents g 

30 npixm=npix ! " : - 

nhxinbands 9/8/g2 , 

deltime=Dtime (timearray) 


Calc histogram min.max of input 

~"wrlte(6,100)infil - 0) 
00 foraatC input file is ..50 ^ iJnage . 

„ote: Tnis min/n-ax call is - reared 

^^^^^SS'SiRix-i- completer 

write ( o» ) 1 

^bSsiS^sSSasssi - - — ' — — 

write (6,*) l 


45 write (6,*) 

write(6,150) GLmin . GL*ax ) 

150 format (/ . 

Do ib«l,nbands 
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75 


write (6, 151) ib, grain (ib) , gmax (ib) 
151 f ormat (5x, 14, 2i9) 

Enddo 

..test: if image mins & maxs are within LDT range of 0-4999: 
Do ib=l,nbands 

If (gmax(ib) .gt. 4999) then 

write (6,*)' Warning — gmax .gt. 4999 STOP 

Call Exit(l) 
Elseif (gmin(tb) .It. 0) then 

write (6,*)' Warning — gmin .It. 0 STOP' 

Call £xit(l) 

Endif 
Enddo 

Note: If gmax gt 4999, then either rescale input data, or redimension 
his arrays in selav3 & selavb, & reset threshold on above tester. 


Create EPS-filtered image 


write (6,*)' [Beginning filtering]' 
If (np .ge. 3 .and. np .le. npm) then 

Call Winbust 042193 (infil, outfil, ibuff, buff n, 

1 npixm, nEx, nlm, np, nl , groin, gmax, pixf orroo , britelo, bnteni, 

2 diagflg) 

Else 

Write (6,*) ' Np out of range stop' 

Endif 


Create EPS residuals image 


25 call Add Image (inf 11, outf il, resf il, f actl, f act2, bias, pixf ormro) 

Treslduals file is unsigned byte w/ mean="bias" ! 

write (6, 200) outf il 
200 formate Ending; EPS file is f ,a50) 

write<6,201) resfil a e n » 
formate EPS residuals file is ',a50) 


201 

30 


deltime=Dtime (timearray) 
write (6, 800) timearray (1) , timearray (2) 
80C format </' User time (sec):', 

x f8.2,'; system time (sec) :' f f 8-2/) 

end 


35 


Ibismnmx: (5-28-91) Calcs. min & max values in a multiband histogram. 

Subroutine Ihismnrax (hist , nbands , gmin, gmax) 
40 Integer hist (5000, 3) , nbands, gmin (3) , gmax (3) 

Do ib=l, nbands 

Do x=l,5000 

If (hist(k,ib) .gt. 0) then 
gmin (ib) =)c -1 
go to 50 

45 Endif 
Enddo 

50 Do x-5000, 1, -1 

If <hist(*,ib) .gt. 0) then 
gmax (ib) ~k -1 
go to 60 


Endif 

Enddo 

60 continue 
Enddo 
return 
end 
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Job: read hdr isl.vf 
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- Read Hdr ISL: (9-25-90) Reads key info from header of an ISL-frat 
\ - image header. Closes image^f ile^upo^conclusion^^^_ = 

5 subroutine Read_Hdr_ISL (filename, npix, nlin, nbd, pxfonn) 

Integer imgnum, ios, npx, nln, nbd, pxf orm, rwmode 
IntegerM opnimg, getdef 
Integer*2 npix, nlin 

•Character*50 filename _ 

i — — *~ 

io ; include ' /local/include/iopacxagef ' 

rwmode =0 

imgnum=0 _ 
ios=Opnimg (imgnum, filename, rwmode, .FALSE. ) 
If (ios -NE. SYSNRM) then „ ^ ^ , 

vrite<6,*) ' Error during 0PNIMG-- Return' 
75 return 
Endif 

ios=Getdef (imgnum, npx, nln, nbd, pxform) 
If (ios .NE. SYSNRM) then 

write (6,*)' Error during Getdef — Return 
. return 

Endif 

Call Cls img( imgnum) 
npix^npx 
nlin=nln 

return 
end 
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70 


8 

• Ihister- (5-28-91) Calcs. multiband histogram of an image using a line 
; buffer only. Intended for use w/ select ave2 .for . 

Subroutine Ihister (infil , npixm, nbx, ibuf f , hist ) 
Character*50 infil 

inclHr SfS^status. gstatus. r status, imgin. npx. nln. nbands. pixf on. 
integer opnimg, getdef , rdline, line 
integer hist (5000, 3) 
integer*2 ibuf f (npixm, nbx) 

include ' /local/include/iopacJcage . inc' 
arrtyp=IDINT2 

-....Initialize the histogram: 
Do ib~l,nbx 

Do k=l, 5000 

hist (k,ib)-0 

Enddo 

Enddo 

. . en the image; 
imgin~0 

ostatus=opnimg (imgin, infil, 0, . false . ) 
If (ostatus .ne- SYSNRM) then 

. write (6, 15) infil cn t ^ , , 

15 formate Ihister: input image ',a50,' not open—abort ) 

stop 

Endif 

gstatus=»getdef (imgin, npx, nln, nbands, pixform) 
If (gstatus .ne. SYSNRM) then 

write (6,20) . . o , ^ , 

20 formate Ihister: Gstatus error on input— abort ) 

stop 
Endif 

: Read the image lines into the histogram: 

Do j=l,nln 

30 rstatus=rdline (imgin, line, -1 , ibuf f , npixm, arrtyp) 

If(rstatus .ne. SYSNRM) then 

50 format V Ihister: error reading line ', i4, ' —abort' ) 

stop 

Endif 

35 Do ib=l, nbands 

Do i=l,npx 

k=ibuff (i f ib) +1 

hist <k,ib]=hist Oc,ib> +1 

Enddo 

Enddo 

Enddo 

40 

Call clsirog (ungin) 

return 

end 
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10 


Winbust 042193: . j 

Uses a barrel-shifting line buffer in multiple bands. 
Implements a "fast-histogram" method. 


?5 


20 


25 


30 


7 infil: Ch*50, name of input image file; not overwritten, 
outfil: Ch*50, name of output filtered image file. 

Int*2' all-band line buffer array; must be dimensioned by 

calling program. m 
Int*2, strip buffer array; must be dimensioned by calling 

program. . 
Int*4, Max. allowed width of image, m pixels per line. 

In?* 4 1 Actual width of image, in pixels per line. 
Read from image header. V««.,- 
Int*4, Max. allowed no. of bands in input image. Inpu^. 
Int*4, Max. allowed height of filter window (in lines) 
& thus max. allowed height of strip buffer. 

Int*4, Actual "width, height of filter window in pixels, 
lines. Input. 

, 3-eleraent vectors. Assumed min & max band GLs 
of input image. Assumed previously computed. Input. 
Note this limits current max. no. of bands to 3. 
pixformorlntM, Datatype of output filtered IPP image file. 
* E.g. 3=uns byte, 4=signed Int*2, 6=signed R 4. 

Input. 

gmn,gmx: IntM. Returned via Common. Set to cur rent „ a1lloe 

-band value of gmin & gmax, respectively. These values MUST 
be set prior to calling Iselavh; Function Iselavh uses 
these values passed via Common.. *--,«.^ ^„h,-« 

nav: lnt*4, no. of samples to average in each filter window. 
Passed thru Common by caller. 

General comments: 


ibuff : 

buffn: 

npixm: 

npix: 

nbx: 
nlm: 

np, nl : 
gmin , gmax : Int * 4 , 


(1) 


(2) 


99% of the complexity of this subroutine is due to 
the implementation of the strip^buffer. The advanta|e of the 
strip buffer is that it allows -convolution"-type window 
operations while keeping at any one time in ™ e ™ or Y °"£y nH - for 
npix*nl*nbands worth of image, instead of npix*nlin*nbands for 

io?rtha? U tSe^s a li?ter is *N0T* recursive, lie. it only uses 
OLD? original GL values in its input weights; . thus output from 
previous window filter do not a^^igtoorw^indwjut^ 


35 


40 


45 


Subroutine Winbust 042193 (infil, outfil, ibuff , buff n, npixm, nbx, 
nlm, np, nl , gmin, gmax, pixf ormo, britelo, britehi, diagf lg) 

Character*50 infil, outfil 

integer imgin, imgout , ostatus, gstatus 

integer nbx, nbands , pixf orm, pixf ormo . noixb 

integer pstatus, ostatus2, r status , wstatus, npx,nln, npixm, npix, npixb 

integer nlin,nlinb 

integer arrtyp, lininc, jtype 

integer ginn, gmx, Iselavh, nav, gmin (3) , gmax (3) 

integer nlm, np, nl , nl2 , np2 , nbinx , gl_cent 

integer opnimg, getdef , putdef , rdiine, vrime 

integer diagf lg 

integer* 4 britelo, britehi, count 

integer*2 ibuff (npixm, nbx) , buff n (npixm, nlm, nbx) 
integer*2 bis (5000) 


50 


55 
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11 

integers els (4096) !Note implicit max. of npix=4096 ! 

data nbinx/5000/ 
common /par/cmn, gmx, nav 

iizx;? 1 '" 2 ssfis-s-s as? t d0 " 

np2-np/2 ! Half-width 

count=0 

"r^en'inputToucpurimages for sequential I/O: 
imgin=0 

ostatus=SYSNRM . ^ > 

ostatus=opmmg (imgxn, mnl,u, .xaxse. ; 
If (ostatus .ne. SYSNKM) then 

?^Itr 7 SUCs^^ w '' a50 <' not w* 0 "'* 

20 Call Exit(l) 

Endif 

SS^tu^S^f {i «gla. npx, nln, nbands, pi«f ox.) 
If (gstatus -ne. SYSNRM) then 

f^t'licatus error on input -abort', 

/j Call Exit(l) 

Endif 

noix-npx ! Actual no. pixels/line in image 

gi^ix W^^-o^'pS.JSSS in filter buffer 
If (npi ^ r iS(6?*r m winbus t _0 4 2193: Npixb > Npixm-Abort' 
Call Exit(l) 

Endif 

ft&SSi&tSiTK. of bands in image 


25 


35 - ;5^5>?°ininbu.t 042193-: «Pix Hlin Type Nbands Typeo N P NaV . 

;r^e^?8^ t Six D n?In^ilcfom.nband S ,pixfor,no.n P ,nav, 


80 


1 
1 


81 
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britelo f britehi r diagflg 
format <21x,2i5, i6, 17 , 16, i4, 15,2x7, 18) 

imgout=0 

pstatus^de f a»™n P x.nln,nbands, P ixforn 1 o, 

?f ^"^tc^t-.f'^^rror for output-abort' 
Call Exit(l) 

Endif 

OStatUs2-SYSNRM ^ a i« 1 

ostatusa-opnimg ^^fjifci^g!'- ' 

If (OS "^te-(6?i)"Sut t irage not open-abort' 
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Endif 


Call Exit(l) 


12 


lininc=nlin/8 IVariable used in writing current line-number 
jtype-lininc ! Ditto (just screen i/o related) 
rstatus=SYSNRM 
wstatus=SYSNRM 


70 

110 


z Load 1st strip into barrel buffer: 

Do 1-1, n!2 +1 !Do per line in 1st buffer strip: 

line=j -1 ! (Damned "zero-starting" C array convention) 

rstatus=rdline (imgin, line, -1, ibuf f , npixm, arrtyp) 
If (rstatus .ne. SYSNRM) then 


75 


vrite(6,110) j , . M , - . -, . 

formate "Winbust_042193" : Read error at line ',i4,' -abort ) 
Call Exit(l) 
Endif 

! Transfer GLs from line buffer to strip buffer: 

! (Someone should be able speed up all this moving around of data) 
!i =horiz pixel index for original image; . 
!i2=horiz pixel index for central (non-mirror) pixels in strip bi 
!i3=horiz pixel index for edge (mirror) pixels in strip buffer. 

j2=j +nl2 
j3=nl2 +2 -j 
Do ib=l,nbands 

20 Do i=l,np2 +1 !Left 

i2=i +np2 
i3=np2 +2 -i 

buf fn <i2, j2 , ibj =ibuf f (i, lb) 
buffn <i3, j2, ib)=ibuff <i, ib) 
buf fn (i2, j3, ib) =ibuf f (i, ib) 
buffn (13, j3 , ib) *ibuf f (i, ib) 

25 Enddo 

Do i=np2 +2,npix -np2 -1 'Central 
i2=i +np2 

buffn (i2, j2, ib)=ibuff (i, ib) 
buffn (i2,j3,ib)=ibuff(i,ib) 

Enddo 

. Do i=npix -np2,npix ! Right 

i2=i i-np2 
i3=2*npix +np2 -i 
buf fn (12, j2, ib) =*ibuf f (i, ib) 
buffn (i3, j2, ib)=ibuff (i,ib) 
buffn(i2, j3,ib)=ibuff (i,ib) 
buffn (i3, j3, ib)=ibuff (i, ib) 

Enddo 


30 


35 


Enddo 

•Write a status display to standard output! 

If <(j -nl2 -1) .ge. jtype) then . 

write (6, 600) j !Write current line number to stdio. 

jtype=jtype +lininc 


Endif 

40 Enddo !End top -j* loop! 


45 


jend=nl !"jend" is barrel index to most recent line 
jc-jend -nl2 !"jc" is barrel index to ■middle" line 
Above sections only serves to fill the strip buffer for the 

first time; it gets more complicated as the buffer contents move 
"down" in the image. 
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w 


30 


the buffer line index is rotated. 

■ T-j-vac rh*> (*last* line +1) of current 
Do j«n!2 +2,nlin +nl2 +1 . ! Indexes tlie C i«t 

J j 0= j - n i2 -1 ! Image index of output lane. 

•Mark which pixels in the current line will be filtered: 
Do i=np2 +l,npixb -np2 

- ""^BfSfcfcS! :S: — — 


Endif 

Enddo 


cls(i)=0 
go to 210 


r-lsfil=l !mark for filtering! ^ t _ z . 

countUount +1 ! Increment to-be-filtered pixel counter! 

75 210 continue 
Enddo 

Do ii>=l , nbands 

gmn»groin ( ib) 
gmx=gmax ( ib ) 

20 i Load window in central P a ^_°:_**ff fL--f £ffl 

! ]? l^VSS? ^^^^^^ band 

1 N ot2°thtt r £xe!r within i/2-windov width of vert image boundaries 
25 j Be ?° w Sii'SSS'as above -12". i.e. strip buffer horiz index. 

._ n ., !Center-of-window index, 1st window! 

iS! ! Output pixel index! 

Call Hist init <nlm,np,nl,np2,buffn,npixn,nbx, 
npixb, ribands, i, jend, ib, nis, nbinx) 

g l_cen t =buffn<i /jC .ib, i^^^fo" . 

" <ClS< i b uf! < ?ioHb^Iselavh(his,nbinx,gl_cent) -Filter 

Else ^ ' 
ibuff (io, ib) =gl — cent 

35 Endif 

! update histogram for later window locations in a line: 

Do i=n?2 +2,npixJb -n P 2 ! Cent er-of -window index! 
io=i -np2 ! Output pixel index! 

40 ~ Call Hist update <nlm, np, nl, np2, buff n, . . , - 

! npixH,nbx,npixb, nbands, ioend,ib, his, nbinx) 

tif passes brightness thresholds test, then filter: 
gl_cent-buffn<?, jC ,ib> SVLii'lS.aS? 

If <C i Du ^ ( l|?ibtilselavh(his,nbinx,gl^^ ! Filter 

Else 
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ibuff <io, ib)«=gl_cent 
Endif 14 


Enddo 


Enddo !End band loop! 

line=jo -1 

wstatus=wrline (iragout, line, -1, ibuff, npixm, arrtyp) 

If (wstatus .ne. SYSNRM) then 

write (6, 410) jo u _,v 

formate *Winbust_042l 93" : write error at line ',i4, -abort') 
Call Exit (1) 

Endif 

♦Rotate barrel-buffer line indices: 

jend=jend +1 

If (jend .gt. nl) jend=l 

jc=jc +1 

If <jc .gt. nl) jc=l 

If (jo -le. nlin -nl2 -1) then 

j2=jend IFor bottom mirror! 

IKead new image line: 
line-j -1 

rstatus^rdline (imgin, line, -1 , ibuff, npixm, arrtyp) 
If (rstatus .ne. SYSNRM) then 
write (6, 430) j 

formate "WinbustJ^lSS" :Error reading line i4, -abort 
Call Exit (1) 

Endif 

Do ib=l, nbands ! Overwrite oldest line in barrel-buffer! 
Do i=l,np2 +1 !Left 
i2=i +np2 
i3-np2 +2 -i 

buffn (i2, jend, ib) =ibuff (i, ib) 
buff n (i3, jend, ib) =ibuff (i, ib) 

Enddo 

Do i=np2 +2,npix -np2 -1 ! Central 

i2=i +np2 

buffn (i2, jend,ib)=ibuff (i,ib) 

Enddo 

Do i=npix -np2,npix ! Right 

i2=i +np2 
13=2*npix +np2 -i 
buffn (i2, jend, ib) =ibuf f (i, ib) 
buffn U3, jend, ib) =ibuf f <i, ib) 

Enddo 
Enddo 

Else 'Virtual mirror of bottom of image: 

j2=jend +1 ! Buffer index of "top" line in buffer! 

If <j2 .gt. nl) j2=l 

j2=j2 -1 ! (12-18-92) 

If (j2 .It. 1) j2=nl '(12-18-92) 
Do ib-1 , nbands 

Do i«l,npixb 

buffn (i, jend, ib) -buffn (i, j2, ib) 

Enddo 
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Enddo 

Endif 

t Write a status display to standard output! 
if <j -ge- jtype) then 

write (6, 600) j , . corop i e ted' 

format UOx, ' Line ',x4, is com? 
jtype=jtype -rlinmc 


Endif 

Enddo iEnd main • j" loop! 


write (6, BOO) count ntered- ' , i8) 

^^•)'°" fi?£5££SEt- coveted,' 

Call clsimg (imgin) • 
Call clsimg (imgout) 
return 
end 
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17 

Hist init: (8-12-91) Used in WinbustJ>42193. Initializes S calcs a histogram 
" of a window. 

Subroutine Hist init (nlm, np,nl ; np2 ; buf fn,npixm, nbx^ 
1 npix,nbands, i, jend, ib,his, nbinx) 

integer *2 buf fh (npixro, nlra, nbx) 
integer*2 his (nbinx) 

integer nlm, np, nl , n P 2 , npixm, npix, nbx, nbands , i , jend, ib, nbinx 
integer gron,gmx,jO 

common /par/gmn, gmx 

...initialize histogram to zero: 
Do k«gmn+l, gmx+1 
his(k)=0 

Enddo 

...Load histogram of window: 
j0= jend 
Do jj-l,nl 

jO-jO +1 

If (10 -gt. nl) jO=l 
Do ii«l,np 

i0~i -np2 -1 +ii 

ia»buffn(iO, j0,ib) +1 

his (ia) -his (ia) +1 

Enddo 

Enddo 

return 
end 
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inteqer*2 buf fn (npixm, nlm. nbx) 
integers his (nbinx) 


, « »l n P 2 npixm,nbx,nbands,i,jend,ib,nbinx 
integer npm f nlm,np,nl, PPf' npiXBl ' 
integer jO, idrop, iadd, npix 


jO-jend 

idrop^i -np2 -1 
iadd=i +np2 


Do 1j=l#nl 

jO=jO +1 


<f Tjo Tgt. ni) 10-i 

ia^buf f n (idrop, jO, ib) 
his(ia)-his(ia) -1 
ia=buffn(iadd, jO,ib) +1 
his(ia)«his(ia) +1 

Enddo 


+1 


return 
end 
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Job: iselavh 042193. vf 
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Iselavh: (6-19-91) Mod'd from selavb to eliminate heavy central 21 
weighting on initial-point in interest of speed. 
Should be valid for all window sizes. 

Note that this function runs more slowly as the image dynamic range 
tgmax -gmin] increases. 
Assumes histogram is precalculated. 

"ac" is the GL of the "center- pixel in the current window. 
References : 

(1) "Noise Cleaning by Iterated Local Averaging", by L. Davis 
« A. Rosenfeld, IEEE Transactions Systems, Man, Cybernetics, 
Vol. SMC-8, No. 9, pp. 705-710, Sept. 1978. 

(2) " Digital image smoothing and the sigma filter", 

by J. S. Lee, CVGIP, V. 24, pp255-269 (1983). Sigma filter, however, 
uses a *fixed* range, in terms of +/- a*sigma, where sigma=est. 
std. dev. of additive noise. 

This function computes a selective-averaging filter output for a 
current location of a rectangular window of 1 band of an image array, 
& returns this output as the standard function return. 

Key Parameters: 

his: Int*2 (nbinx) ; Histogram of GLs in current window; 

Assumes that GL value corresponding to a histo bin is 

offset by 1, e.g. his(l) corresponds to GL=0. 
nbinx: Int*4; Max. no. of bins in histogram, 
ac: Int*4; spatially-central GL in current window. 

gran,gmx:lnt*4; min & max GLs in entire image for present band. 

Passed from caller via Common, 
nav: Int*4; desired no. of pixels to average in the 

current window. 

Output "iselavh" will be a weighted sum of histogram bin occupancies; 
this weighted sum is num/den. 


Integer Function Iselavh (his, nbinx, ac) 

real* 4 val, vail, val2 

integer*2 his (nbinx) 
integer* 2 cnt 

integer gmin, groax, nav, a c, nbinx 
Integer ar, al,denr, denl 
INTEGER* 4 NUM,DEN 

common /par /gmin, gmax, nav 

common /par2/num, den,ar, al f cnt, vail, val2,val 


. ! average the nav-nearest grey levels . . . 
den-his (ac +1) 

if (den .ge. nav) then • !If have >den local vals =*ac,then no-change! 
iselavh=ac 
return 

endif 

num«ac*his (ac +1) 
do k-1,4999 

ar-ac +k ! Current GL to "right" of center GL 

al-ac -k !- - "left" of " 

if (ar .le. gmax) then 

num=num +ar*his(ar +1) 
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den=den +his (ar +1) 

endif 

if (al .ge. gmin) then 

num=num +al*his(al +1) 
den =den +his (al +1) 


22 


endif 
if (den 
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else if 


enddo 


Endif 

write < 6, *) 
write (6,*) 
write (6,*) ' 
write (6,*) ' 
write (6, *) ' 
Call Exit (1) 


. gt. nav) then 
cnt=den -nav 
denr=0 
denl=0 

If (ar .le. gmax) denr=his (ar +1) 

If (al .ge. gmin) denl=his (al +1) 

val2-Float (denr +denl) 

vall=Float (ar*denr) +Float (al*denl) 

val=vall/val2 

iselavh=Nint ( (Float (num) -cnt*val) /Float (nav) ) 
return 

(den .eq. nav) then 

iselavh=Nint (float (num) /float (den) ) 

return 


Unexpected termination in iselavh! 
Local Iselavh values at Exit:' 
Central GL value (ac) :',ac 
Groin, Gmax: gmin, gmax 


Nav: 


, nav 


end 
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Add_Image: (6-26-91) Adds 2 1-3 band images w/ additive 
constants* 

Assumes the internal image representation is Real*4 
Max. # pixels/line=4096. 

Subroutine Add Image (infill, inf il2, outf il, factl, fact2, const, pixforrao) 
Character*50 inf ill, inf il2, outf il 
Real*4 rbuf f 1 (4096, 3) , rbuf f 2 (4096, 3) 
Real* 4 factl, fact2, const 

integer Imginl, imgin2, im^out, nimx, ostatus, gstatus 
integer nbands, nbandsc, pixf orm, pixf onnc, pixf ormo, band 
integer pstatus, ostatus2, rstatus, wstatus,npx, nln, npx2, nln2 
integer opnimg, getdef , putdef , rdline, wrline 
integer arrtyp, llninc, jtype, imax 

data imax/ 4 096/ 

include ' /local/include/iopaclcagef ' 
arrtyp* 1 IDREA4 

...Open input color image: 
imginl=0 

- ostatus=opnimg (imginl , inf ill, 0, .false. ) 
If (ostatus .ne. SYSNRM) then 
write (6, 10) infill 
10 formate In color image ' r a50,' not open-abort') 

Call Exit(l) 

Endif 

gstatus»getdef (imginl , npx, nln, nbands, pixf orm) 
If (gstatus .ne. SYSNRM) then 
write(6,20> 

20 formate Gstatus error on input — abort') 

Call ExitU) 

Endif 

If (npx .gt. imax) then 

write (6,*) ' Npx .ge. imax — abort' 

write (6,*)' Npx=',npx 

write (6,*)' lmax= ',iniax 

Call Exit(l) 
Endif 
npix=npx 
nlin=nln 
write (6, 45) 

45 format (/' "Add Image": Npix Nlin Type Kbands Typeo' ) 

write (T, 50) npix, nl in, pixf orm, nbands, pixf ormo 
50 format (13x,2i5,i6,iB,i8) 

...Open 2nd input class map: 
imgin2*0 

ostatus-opnimg (imgin2, inf il2, 0, .false. ) 
If (ostatus .ne. SYSNRM) then 

write (6, 60) infil2 
60 formate In 2nd image ',a50,' not open-abort') 

Call Exit(l) 

Endif 

gstatus*getdef (imgin2,npx2,nln2, nbandsc, pixf ormc) 
If (gstatus .ne. SYSNRM) then 
write(6,65) 

65 formate Gstatus error on 2nd input — abort') 

Call Exit(l) 
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3 


"Ife i . n !^r^S 2 aiSIns?^L 5Ku.l-C.ll Exit a," 

Call Exit (1) 
Endif 

pstatus^lSefamgout.npx.nln.nbands.pixformo) 

Call Exit(l) 
w?itus2=opnimg (imgout^out "l^i^fj f 5 ' ' 
„ (o^atus^.ne^SYSNPM^then^ ^ 

Call Exit (1) 

Endif 

write (6, * > ' Factl, Fcat2. Const=' , factl, fact2. const 

.Read each line of input image: 
lininc«nlin/4 
jtype=lininc 
Do j-l,nlin 

line=j-l 

i for all bands! 

If (r status .ne. SYSNRM) then 

r o ^trtoor reading image line -abort') 
Call Exit(l) 
Endif 

«ss«s^£^^T s "" 2 ' li, - i ' arrtw) 

i for all bands! 

If (rstatus .ne. SYSNRM) then 

?ormIt 7 ( ;8 |rror reading image line -abort') 
Call Exit(l) 
Endif 

Do band=l , nbands 
Enddo 

Clio data to output range: 
"if (Pixfonno .eq. IDBYTE) then 

Do i=1 ' n S uffl< i. b and)-toaxl(rbuffl(i.band).0^ 

rbuff 1 (i.band) -Aminl (rbuf fl (1. band) ,255 .) 

Blse^tf" (pixfonno .eq. IDINT1) then 

Do i=l. npix band) , a> , 128 ) 

rbuff 1 (i.band) -Aminl (rbuf fl (i.band) .127.) 

Else^lf 0 (pixfonno .eq.. IDINT2) then 

Do ^^^^(i.bandj-Anintlrbufflj^band)) 

rbuf f 1 (iiband) -Amaxl (rbuf f 1 (i.band) , -32768 . ) 


5 
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100 


w 


rbuf f 1 <i, band) =Aminl (rbuf f 1 (i. band) , 32767 . ) 26 

Enddo 

Else if (pixformo .eq. IDREA4) Chen 
Else 

format r 1 pixfonno=' ,13, ' not implemented-ar>ort ' ) 
Call Exit(l) 
Endif 

Enddo !End band loop! 


c Write output line of data: 

wstatus=Wrline (imgout , line, -1, rbuf f 1, imax, arrtyp) 
If (wstatus -ne. SYSNRM) then 

write (6,*)' Error in wstatus-abort ' 
. . Call Exit (1) 
Endif 

If (j -ge. jtype) then 

write (6, 200) j . 
'5 200 format <10x, ' Line ',14,' is completed ) 

jtype=jtype +lininc 

Endif 

Enddo !End line loop! 

c Close images 

Call cisimg (imginl) 
20 Call cisimg <imgin2) 

Cal 1 els img ( imgout ) 
return 
end 
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l 

Bust Class Creates ** two-state ** class map from 3-band residuals image. 
"Optionally dilates classmap before output. T «. ^ * i 

Calls: Irgb_load_isl, Buster_Class_Create, Dilate_Class, Ioutput_isl. 

Image arrays: . . 

res :3-band dirty residuals image (input) 

els :l-band classmap (output) 

residuals image & class map are assumed lnteger*2 values (i.e., not rea 
Note demonstration max image sire of 1870 x 1374. 

Does not assume that residuals image is biased by a predetermined 
constant; this info is implicit in the GL-res thresholds entered by the 


Calls: 


Irgb load isl . 

~ (calls image-format-specific 10 routines) 
Buster_Class_Create 
Dilate Class 

Ioutpu^ isl , * 

"(calls image-format-specific io routines) 


Program Bust_ Class 

Integer*2 res (3, 1870, 1374) 

Integer*2 els (1870, 1374) , cls2 (1870, 1374 ) 

integer*2 npix,nlin 

Integer npixm, nlinm,band, nbands, nbc, nbx, cnt , np, nl 
integer pixf ormi, pixf ormo 


20 


integer cnttot 
integer ixO, iyO, nx, ny 

integer dgl, glmin, glmax,gltmin, gltmax, dilf Ig, byte 

real* 4 deltime, Dtime, timearray (2) 

character*50 inf il, outf il, resf il, clsf il 
character* 1 iname (50) , iname2 (50) 
character*50 commSO 

equivalence (comm50, iname2 (1) ) 

data dgl/255/ 'Defines "dirt" class! 

data byte/3/ 

data nbx, npixm, nlinm/3, 1870, 1374/ 

data pixformo/3/ ! Classmap of type "byte 

i — — — - ~ 

write (6,*)' IMax. image size: 1870 x 1374)' 
write (6,*)' ' 

format renter in biased residuals image filename [ch50]:',$> 
read (5,' (a50) ' ) resfil 


write (6, 30) 

30 fonnatC Enter out classmap image name: ,9) 

read(5, # (a50) ' ) clsfil 

write(6,50) w w _ , e . 

50 formate Enter min, max residuals GL thresholds: ,5i 

read (5,*) gltmin, gltmax 

write (6, 60) 
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2 

60 " format <' Do you wish prefilter dilation? [yes="l"J : ' , $) 

read (5,*) dilflg 

c - 

5 write (6,*) ' Class threshold used is ' , dgl 

deltime=Dtime (timearray) 

!Read biased color-residuals image: 
10 Call Irgb_load isl (resfil , res, nbx, ribands, npixro, nlinm, 

1 npix, nTin,pixformi, 0) 

write (6, 150) resfil 
150 formate Input residuals file is ',a50) 

t — 

;s ! Compute class map: 

Call Buster_Class_Create (res, els, nbx, nbands, 
1 npixrn, nlinm, npix, nlin, cltmin, gltinax, dgl) 

i — — 

!Optionally, dilate classmap by one pixel: 

20 If (dilflg .eq. 1) then 

write (6,*)' [Dilating class map]' 

Call Dilate_Class (els, cls2, 
1 npixm, nlinm, npix, nlin, dgl) 

25 Call loutput_isl (clsfil, cls2, npixm, nlinm, 

1 npix, nlin, byte, inane) 

write(6,200) clsfil 
200 formate Out dilated class file is ',a50) 

30 Else 'write undilated classmap! 

Call Ioutput_isl (clsfil, els, npixm, nlinm, 
1 npix, nlin, byte, iname) 

write (6, 300) clsfil 
300 formate Out undilated class file is ',aS0) 
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Endif 

deltirae=Dtime (timearray) 
write (6, 800) timearray ( 1) , timearray (2) 
800 format (/' Bust_Class:: user time (sec):', 
1 f8.2,'; system time (sec) : ' , f 8 . 2/) 


end 
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c irgb_load_isl: Loads up to a 3-band image into a Int*2 array. 

Subroutine Irgb_load_isl (bname,b,nbx,nb, npixm,nlinro, npix, nlin, pixform, 
1 lunprt) 
5 Characterise bname, infil „ ilftfl , 

Integer*2 b (nbx, npixin,nlinm) , npix, nlin, ibuff (4096, J) 

integer imgin, gstatus, ostatus 

integer nbands, pix form, band, nbx, nb,npixra, nlinm 

integer rstatus,npx,nln 

integer opnimg, getdef , rdline 

integer arrtyp, lininc, jtype,pixdim 

10 

Data pixdim/4096/ 

include ' / local / include/iopackagef ' 

arrtyp=IDINT2 

c Open input image: 

;s imgin~0 

inf il=bname 

ostatus-opnirog (imgin, inf il, 0, . false . ) 
If (ostatus .ne. SYSNRM) then 

write (6, 10) infil 
10 formate In image ',a50,' not open-abort' ) 

stop 

20 Endif 

gstatus^getdef (imgin, npx, nln, nbands, pixform) 
If (gstatus .ne. SYSNRM) then 

write (6,20) w _ t 

20 formate Gstatus error on input — abort ) 

stop 

25 Endif ^ % _ 

If (nbands .gt. 3 -or. nbands .gt. nbx) then 

write (6,*)' Input no. bands too large — abort 
stop 

Endif . ^ 

If (npx .gt. npixm .or. npx .gt. pixdim) then 

write (6,*)' -IRGB_load_isl n : too many pixels— abort 
30 stop 
Endif 

npix=npx 
nlin=nln 
nb«nbands 

write (6, 45) npix, nlin, nbands, pixform . M . . 9 

35 45 format (2x, ' "Irgb load_isl":: Npix:',i5,'; Nlm:',i5,'; Nb: ,i2, 
1 '; Type:'7i2) 

If (lunprt .gt. 0) then t 

write (lunprt, 45) npix, nlin, nbands, pixform 

.Endif 

c Read each line of input image: 

40 lininc«nlin/4 
jtype-Xininc 
Do j-l,nlin 

c Read input lines & load 3-d array: 

line=j-l , Jt 

rstatus=rdline (imgin, line, -1 , ibuf f , pixdim, arrtyp) 
i f 0r a n bands! 

45 If (rstatus .ne. SYSNRM) then 

write (7, 80) j 

80 formate Error reading line -abort') 

Stop 
Endif 
Do 1-1, npix 

50 
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200 


Do band=l,nbands 

b (band, i, j) =ibuf f (i,band) 

Enddo 
Enddo 

If <j .ge. jtype) then 
write<6,200) j 

format UOx*, ' Line ',14,' has been read') 
jtype= jtype +lininc 

Endif 

Enddo !End line loop! 


10 


.Close image: 

Call clsimg (imgin) 

return 

end 
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Job: buster class create_dugas.vf 
Date: Wed Apr 21 17:46:32 1993 
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c Buster aa»sJCre.te S Cre«te8 .^rt^^te^^s^^rom 

c the biased color residuals of a f^tSmootn typ 9. 

% SStrS^lS^rSi?lS&«?5S2« d *.t« of presence 

c of defects. 

Subroutine Buster Class_Create (res, els, nbx, nbands, npixm, nlinm, 
1 npix, nlinTgltmin, gltmax, dgl ) 

Integer *2 res (nbx, npixm, nlinm) , els (npixm, nlirun) 
Integer*2 npix,nlin 

integer nbx, nbands, npixm, nlinm, gltmin, gltmax, dgl, glmin, glmax 


20 


25 


10 
11 


formlM ;1 Bust Class_Create: Gltmin Gltmax Dgl') 
write (6, 11) gTtmin, gltmax, dgl 
format (19x, 2i7, 15/) 

Do j«l,nlin 

Do i-l,npix 

glmin*res (1, i, 1) 


Enddo 

return 

end 


Enddo 


g Imax^g lmin 
Do ib«2, nbands 

glmin«MinO (glmin, res (ib, i,J) 
glmax»MaxO (glmax, res (ib, i, j) ) 

?? d tglmin .ge. gltmin .and. glmax .ge. gltmax) then 
cls(i, j)-=dgl 

Else 

cls(i, j)~0 

Endif 
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Job: dilate class dugas.vf 
Date: Wed Apr 2 1 lT:47:31 1993 
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c Dilate Class: Dilates by 1 ? ixe i * s? 1 ? f *%^ a " ^window, 
c Byte class-map. Dilation kernel is a 3x3 square 


w 


15 


30 


ctar-class which is dilated; if current pixel is undilated then 
old classvalue is copied as new classvalue. 

A1S ° ^'S^SiSfl^if^l^ still covering image. 
^^"^atrcHsTin npix, nlin, ctar) 

Integer *2 clsi <npixro, nlinxn) , clso (npixm, nlinm) ,npix, nlin 

Integer npixm, nlinm, ctar, im, ip, jm, jp 

i — • — 

write (6,10) ctar - . 

10 format </' Dilate_Class : target class xs .1^/) 

Do 1=1, nlin 

jm=j-l 

Im-MaxO (jm, 1) 
jP«j + l 

jp«MinO ( jp, nlin) 
Do i«l,npix 

im=i-l 

im=MaxO (im, 1) 
ip=i+l 

ip=MinO (ip,npix) 
Do jp 
Do ii=im r ip 

■If any neighbor^ctar, then 

! set center=ctar as well: 

If <clsi<ii.jj> -eq. ctar) then 
clso<i,1)*ctar 
go to 100 

Endif 

Enddo 

Enddo ^ ^ 

!No change, so copy to output: 
clso(i, j)«clsi<i, j) 
10 q continue 


Enddo 


Enddo 


35 


return 
end 
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SUBROUTINE: IOUl^rUT isl (10-9-90) 
PROGRAMMER: ELIZABETH J. DONALDSON /R- Gray 


nteger*2 - DESTINATION ARRAY ■ 
INTEGER* 4 r MAX PIXEL DIMENSIONS FOR 
INTEGER* 4 - MAX LINE DIMENSION FOR F 2 
INTEGER* 4 - NUMBER PIXELS IN INFJJE 
INTEGER* 4 - NUMBER LINES IN INFILE 
rxArun „ - w« it w. - INTEGER* 4 - INFILE DATATYPE 
INAME50 - OUTPUT - CHARACTER - (UNUSED) 


DATE: 

9/20/90 

OFILE 

- INPUT - 

FOUT 

- INPUT - 

NPIXM 

- INPUT • 

NLINM 

- INPUT ■ 

NPIX 

- OUTPUT 

NLIM 

- OUTPUT 


"subroutine Ioutput_isl*<of ile, f out , npixro, nlinm, npix, nlin, pixf orm, 
inameSO) 

include ' /local/include/iopackage . inc' 

character of ile*50, inameSO (50) 

Integer*2 fout (npixm, nlinm) 

integer*2 npix,nlin 

iSSSS.* f£S£: SSS: *** <« .«*!-. mm- 

integer * 4 ostatus pstatus , wstatus , npx^nln 

i — ~ 

vrite<6,5) npix, nlin, pixf orm , , . , - 4 

formate I output isl: image header data: np = n± f^' 

' pixf orm = ' , i3T 

arrtyp - IDINT2 
imgout=0 

DEFINE CHARACTERISTICS OF OUTPUT IMAGE 
npx=npix 

nln=nlin . 
©status = putdef (imgout r npx, nln, 1, pixf orm) 
if ( gstatus .ne. SYSNRM) go to 1000 

OPEN OUTPUT IMAGE FOR RANDOM ACCESS- 
ostatus *» opnimg (imgout , of ile, 1, .true . ) 
if ( ostatus ,ne. SYSNRM) then 

write(6,10)ofile , 

formatC Output image' , a50, not open. ) 

write (6, *) ' Exiting — ostatus= ostatus 

go to 1000 

endif 

if( pixform.eq.IDBYTE) then 
do 100 j » l,nlin 

do 200 i - l,npix 

if (fout(i,j)-lt.0) fout(i,j) = 0 
if (fout <i, j) -gt.255) fout(i,j) =255 
continue 
continue 
elseif( pixf orm. eq.IDINTl) then 
do 110 j - l,nlin 

do 210 i - l,npix 

if (foutti, j) .It. -128) fout(i,j) « -128 
if (fout (i,j) -gt. 127) fout U,3) '127 
continue 
continue 
elseif( pixform.eq.IDINT2) then 
do 120 j - l,nlin 

do 220 i • l,npix „. cfl 
if ( fout U.1) .It. -32768) , fou t« A '3>r"3"fi« 
fout <i, j) .gt. 32767) fout <i, j) =32767 
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220 continue 
120 continue 

elseif ( pixforro.eq.IDREA4) then 
GO TO 900 

else 

write (6, 30)pixform 
30 formate IOUTPUT_ISL: Pixform ',i3, ' not implemented.') 

endif 

900 wstatus » wband<imgout,0,fout,npixm,nlinm,arrtyp) 
if ( rstatus -ne. SYSNRM) go to 1000 

1000 continue 

call clsirag (imgout) 

return 

end 


Appendix C — Docket 66,719 
© 1993 Eastman Kodak Company, Rochester, N-Y- 14650-2201 
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Buster Erode am; (11-11-92) Uses dynamic memory allocation on the Sun. 
Calls Hedirt5x5 byte2 nopad; no 2-pixel "dead zone 
around image; all pixels in image are potentially cleaned. 

Intended for use on a Sun workstation running SunOS (e.g. a 
Sparcstation 10) . 

5x5 erosion window is currently hard-wired in. Minimum clean fraction 
of erosion window is nominally .24. 

Cleans RGB image by grayscale erosion 
using an external 2-state class map. 

Calls: 

Irgb load isl 

~~ (Includes image-f ormat-specific calls) 
Binput Isl 

~ (includes image-f ormat-specif ic calls) 
Dedirt5x5_byte2_nopad 

Irgb write isl , v 

~ (includes image-f ormat-specif ic calls) . 

Image arrays : - . J . , ^ 

iml -.Original (dirty) 3-band image (input) 

els :l-band classmap (input) 

im2 : Cleaned 3-band image (output) 

Warning: Once erosion iteration begins, iml & im2 contents 
are *partially*-cleaned image; original is not retained. 

Removes by interpolation all pixels whose class-map GL is .ge. 
a threshold value. Works on a multiband image (up to 3 bands). 
Assumes iterative in/of recursive filtering. 

Uses "adjustable variable-input window size, allows only filtering 

a rectangular piece of the image. 
Image & class map are assumed Integer*2 values (i.e., not ^eal)^ 


Program Buster_Erode_Am 

Pointer <piml,iml), (pim2,im2), (pels, els), (pcls2,cls2) 

Integer*2 iml (1, 1, 1) , im2 (1, 1, 1) 
Byte els (1,-1), cls2 (1,1) 

Integer*2 iml (3, 1828, 1332) , im2 (3, 1828, 1332) 
Byte els (1828, 1332) , cls2 (1828, 1332) 


integer *2 npix,nlin 

Integer npixm,nlinm, band, nbands,nbc, nbx,cnt 
integer pixf ormi,pixformo,pixf ormc 

integer iml_bytes, im2_bytes, cls_bytes, cls2_bytes, tot_bytes 
40 integer itmax, it , cnttot 

integer inflg 
integer dgl,byte 
Integer* 4 iystart, iyend 

real* 4 del time, Dtime, timearray (2) , f min 

45 character*50 inf 11, outf il,clsfil 

data itmax/407 !Max. no. of erosion iterations! 

data dgl/255/ 'Defines "dirf class! 

data dlagflg/0/ 
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data byte/3/ ^ ^ , 

data nbx,npixro,nlinm/3, 1828, 1332/ 


write (6^*)' Buster_Erode_Aro : Dynamically allocated memory' 

write<6,10) , , e» 

!0 format ( /' Enter in -dirty" file name[cn5QJ. ,5) 

read (5,' (a50) ' ) inf il 

Call Read hdr isl (inf il, npix f nlin,nbands, pixformi) 
write(6,15)npTx,nlin,pixformi,nbands 
15 formate npix=',i4,'; nlin-',i4,'; type=',i2,'; bands=',i2) 

npixm=npix 
nlinm=nlin 
nbx=nbands 

iml bytes*2*nbx*npixm*nlinm 

im2~"bytes=iml_bytes 

cls~bytes=npixm*nlinm 

cls7 bytes=cls bytes , , _ . _^ 

tot Bytes=iml Bytes +im2_bytes +cls_bytes +cls2 bytes 
write (6,*)' "Memory allocation requested (bytes): ',tot_bytes 

write (6, 20) 

20 formate Enter in class filename:',?) 

read<5,' (a50) ' ) clsfil 

write (€,30) am , WCA , , e . 

30 format ( ' Enter out "cleaned" filename [ch50] : ,5) 

read (5,' (a50) ' ) outf il 

40 formate Writer reqd cleaning fraction -fmin" [Real]:',$) 

read (5,*) fmin ^ _ A _ 

If (fmin .le. 0. .or. fmin -gt. 0.5) then 

write (6,*)' Fmin out of range — abort 

Call Exit(l) 

Endif 

c. ....... - 

write<6,*) ' Class threshold used is ',dgl 

■Read color image to be filtered: 

write (6, 150) inf il 
150 format (/' Input (dirty) file is ',a50) 

deltime*Dtime (timearray) 

piml«malloc(iml bytes) 

If (piml .eq. 0T then , ' J ^ 

write(6,*)' Malloc for -iml" failed— abort' 
Call Exit(l) 

Endif 

Call Irgb_load isl (inf il, iml, nbx,nbands, npixm, nlinm, 
1 npix, nXin, pixformi, 0) 

pixformo-pixformi 'Out datatype set equal to in type! 

write(6,160) clsfil 
160 formate Input class file is ',a50/) 
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pcls=roalloc (els bytes) 
If (pels .eq. OT then 

write <€,*)' Malloc for "els" failed— abort' 

Call Exit (1) 

Endif 

Call Binput isl (clsf il, els, npixm, nlinm, npix, nlin, pixformc, 
1, ru5c) 

pim2«malloc (im2 bytes) 
If (pim2 .eq. oy then 

write (6,*)' Malloc for "im2" failed—abort' 

Call Exit (1) 

Endif 

pcls2 ss malloc (cls2_bytes) 
If (pcls2 .eq. 0) then 

write (6,*)' Malloc for "cls2" failed—abort' 

Call ExitU) 

Endif 

i 

inflg=l 'Start by assuming "els" is input class map! 

! Clean by erosion iteratively until no more pixels are replaced: 

write (6,*)' [Beginning cleaning]' 
cnttot=0 

iystart=l IRequired initializations 

iyend=nlin 

Do it«l,itniax 

cnt=0 

If (inflg -eq. 1) then 

Call DedirtSxS byte2_nopad (iml, els, im2, 

cls2, npixm, nlinm, npix, nlin, nbx, nbands, 

dgl, cnt, fmin) 
Call DedirtSxS byte2_nopad2 (iml, els, im2, 

1 cls7, npixm, nlinm, npix, nlin, nbx, nbands, 

2 dgl, cnt, fmin, iystart, iyend) 
inflg=2 

Else if (inflg .eq. 2) then 

Call Dedirt5x5_byte2_nopad (im2, cls2, iml, 

els, npixm, nlinm, npix, nlin, nbx, nbands, 

dgl, cnt, fmin) 
Call Dedirt5x5_byte2_nopad2 <im2, cls2, iml, 

1 els, npixm, nlinm, npix, nlin, nbx, nbands, 

2 dgl, cnt, fmin, iystart, iyend) 
inflg=l 

Else 

write (6,*)' Inflg out of range — abort' 
Call ExitU) 

Endif 

cnttot=cnttot +cnt 
write (6,200) it , cnt 

formate After it=',i3,' « pels cleaned = ',i6) 

If (cnt .eq. 0) then !No more dirt-class pels to be cleaned! 
go to S00 'Dreaded "go-to" statement! 

Endif 
Enddo 

write (6,*)' Warning-Cleaning loop exhausted after it=',itmax 
write (6,*)' Aborting due to loop exhaustion' 
Call Exit(l) 
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continue 

write (6 r 550) it , 
formats Cleaning ended after iteration^' , i3> 
write (6, 552) cnttot . 
foraatC Total f of points f iltered=' , 16) 

! Release classmap memory: 
Call Free (pels) 
Call Free(pcls2) 

If (inflg -eq. 1) then , 

Call Free(pim2) ! release unr.eeded image memory 
call Irgb_write_isl<outf il, iml, nbx,nbands, npixm r nlinm 
npix, nlin, pixf ormo, 0) 

Call Free(piml) ! released unneeded image memory 
call Irgb write_isl (outf il, im2, nbx, nbands, npixm, nlinm 
npix, nlin, pixf ormo, 0) 


Else 


Endif 


deltime=Dtime (timearray) 
write (6, BOO) timearray (1) , timearray (2 ) 
format (/' Buster Erode_Am: : user time (sec): ( 

f8„2,'; system time (sec) : ' , f 8 .2/ ) 

write (6, 850) outfil . 
formate Cleaned output image file is ',a50) 

end 
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c Binput_isl: (10-7-92) Reads one band of an isl-format image into 
c a single byte-type array. Input is clipped to range 0-255. 


Subroutine Binput isl (inf il, im, npixm, nlinm, npix, nlin, pixf onn, 
1 band, nbands) 

5 

Character*50 infil 
Character cval 

Byte im (npixm, nlinm) , bval 

Integer*2 ibuff (4096) , npix, nlin 

10 

integer imgin, ostatus, gstatus 
integer opnimg,getdef , rdline 
integer nbands, pixf orm, band, iband, ival 
integer rstatus, npx, nln, npixm, nlinm 
integer arrtyp, lininc, jtype 

'5 Equivalence (cval, bval) 

i — 

c include ' /local/include/iopackagef ' 

! Above comment -out must be removed after f77cvt conversion!!!! 

arrtyp=IDINT2 

20 c.-.Open input image: 

imgin^O 

ostatus=opnimg (imgin, infil , 0, . false . ) 
If (ostatus .ne. SYSNRM) then 
write (6, 10) infil 
10 formate In image ',a50,' not open-abort') 

Call Exit(l) 

25 Endif 

gstatus^getdef (imgin, npx, nln, nbands, pixf orm) 
If (gstatus .ne. SYSNRM) then 
write (6,20) 

20 formate Gstatus error on input — abort') 

Call Exit(l) 

30 Endif 

npix^npx 
nlin=nln 
write (6, 45) 

45 format (/' "Binput isl": Npix Nlin Type Nbands Band') 
write(6,5TX) npix, nlin, pixf orm, nbands, band 
35 50 format (14x,2i5,i6,iB,i5) 

c — Read each line of input image: 
lininc=nlin/4 
jtype «=lininc 
rstatus«0 

40 Do j-=l,nlin 

c Read input lines & calc. sum: 

line»j-l 
iband-band -1 

rstatus=rdline (imgin, line, iband, ibuff, 4096, arrtyp) 
If (rstatus .ne. SYSNRM) then 
45 write (7, 80) j 

80 formate Error reading line -abort') 

Call Exit(l) 
Endif 

Do i»l,npix 
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Enddo 


Enddo 


ival-MaxO (ibuf f (i) , 0) 
ival^MinO (ival, 255) 
cval=Char (ival) 
ira(i, j)=bval 

!End line loop! 


.Close images 

Call clsimg (imgin) 

return 

end 
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c Dedirt5x5_byte2_nopaa: 9 

c Eliminates "dead zone" 

c of 2 pixels around image edges (i.e. cleaning now takes place at all 

c pixels of the image) - Only need one clean pixel in a 5x5 neighborhood t 

c clean a dirty border pixel (i.e. nmin_edge=l) . 

5 

c Output non-dirt or cleaned classes in classmap have a value of "0". 

c Note that within one call of Dedirt, cleaned Ghs are *not* used 

c to clean their neighbors; this will, however, occur in 

c subsequent iterations (further calls to dedirt) . 

c Ad hoc "dirt" remover. Assumes that dirt is detectable as 

70 c a simple grey-level threshold in the input classmap. 

c Filters all pixels *above* a certain 

c threshold grey level ("gth") . 

c Note that a 1-pixel dilation of the classmap is suggested prior 

c to first using this subroutine, tc avoid propagation of dirt 

c "ghosts". 

7 c c 

c Image arrays: (Max dims :npixm, niinm; Actual dims: npix,nlin) 

c iml: Input (dirty) 3-band image 

c els: Input dirt 1-band classmap 

c cls2: Output dirt 1-band classmap 

c im2: Output (cleaned) 3-band image 

c 

c Input: np,nl :Size of filter window. 

20 c dgl : Threshold classmap GL value req'd for cleaning. 

c fmin : Fraction of window needed to be "clean" to be filtered 

c (l=<nmin<np* (nl-1) /2 ). 

c Output :cnt :No. of pixels which were filtered. 

Subroutine DedirtSxS Byte2_Nopad (iml, els, im2, cls2, npixm, nlinm, npix, nlin 
25 l dgl, cnt, fminT 

Integer* 2 iml (nbx, npixm, nlinm) , im2 (nbx, npixm, nlinm) 

Integer*2 npix f nlin 

Integer npixm, nlinm, np,nl, np2,nl2, cnt, dent, 10, ix, jO, jx, dgl 
Integer nbx, nb, num_neigh, nmin, delx, dely 
30 Integer nmin_edge, nmin_cent 

Byte els (npixm, nlinm) , cls2 (npixm, nlinm) ,bdgl 

Real*4 sum (3) ,wtsum, wt (3, 3) , fmin 
■fmin=.24 implies nmin=6 for 5x5 window. 

35 Character cdgl 

Data np,nl/5,5/ 'Erosion window size! 

Data vt(l,l),vtC2,l),wt(3,l>/l.,l.,0.5/ !Wts specific t: 

Data wt <l,2),wt (2,2) ,wt (3,2) /l., 0.707, 0.4472/ 
Data wt (1,3) ,wt (2,3), wt (3,3) /0. 5, 0.4472,0.3536/ 
Data nmin_edge/0/ 
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Equivalence (cdgl,bdgl) 


cnt=«0 
dcnt«0 
nl2~nl/2 
np2«np/2 

nmin_cent=Nint (fmin*np*nl) 
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cdgl=Char (dgl) 10 
write (6,10) 

format (' DedirtSxS Byte2_Nopad: Npix Nlin Np Nl Fmin Nmin Dgl') 
write (6, 11) npix, nTin, np, nl, fmin, nmin cent, dgl 
format <23x, 2i5, 2i4, f 8.3, i5,i6/) 


! Slide window across image: 
Do j-l,nlin 
Do i-l # npix 


!i,j are coordinates of center of window! 
•Test local value of classmap: 

If (cls(i,j) .eq. bdgl) then 'Filter the point 
dcnt=dcnt +1 
j0=j-nl2 
jx= j+nl2 
jO=MaxO(jO,l) 
x=MinO < jx,nlin) 
0«=i-np2 
ix*=i+np2 
i0=MaxOU0, 1) 
ix=MinO (ix, npix) 

!iO f ix,jO,jx are coordinates of extrema of 
! current window 
Do ib=l,nb 

sum(ib) =0 . 

Enddo 
wtsum=0 . 
nura_neigh=0 

!Calc. over local "dirty" window: 
Do j2=j0, jx 

Do i2=i0,ix 

!If neighbor is *not* dirt, then: 
If (cls(i2,j2) .ne. bdgl) then 

!Calc. a distance-weight for that GL: 
!dis= (Float (i -i2))**2 + (Float (j -j2))**2 
!dis-Sqrt (dis) 
!wt~l./dis 
! wtsum=wtsuro +wt 


delx«Int(Abs(i -i2) +1.) 
dely=Int<Abs(j -j2) +1.) 
wtsum=wtsum +wt (delx, dely) 

Do ib«l, nb 

sum(ib)«sum(ib) +wt (delx, dely) *iml (ib, 

Enddo 

! Count no. of non-dirt window pixels: 
num neigh=num neigh +1 
Endif "~ ~ 

Enddo 
Enddo 

JCalc. local value of "nmin": 

If (j .le. nl2 .or. j .ge. (nlin -n!2 +1)) then 
nmin«nmin_edge .'Top or bottom border! 

Else if (i .le. np2 -or. i -ge. (npix -np2 +1)) then 
nmin=nmin* edge !Left or right border! 

Else ~ 

nmin«nmin_cent !Not in a border region! 
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ndif 11 

!If have min. no. of non-dirt neighbors, then 
! replace dirty GL w/ wted ave of nondirty GLs: 
If (nun neigh .gt. nmin) then 'Replace the center GL: 
5 Do i5=l , nb 

im2 (ib, i, j) =Nint (sum <ib) /wtsum) 
Enddo , J . 

cls2(i,j)=0 IDon't filter point again! 
cnt»cnt +1 

Else !Not enough non-dirt pixels within window! 
70 !Keep the old (dirty) GL for now: 

Do ib=l, nb 

im2 Ub, i, j> =iml (ib, i, j) 

Enddo 

cls2(i, j)=cls(i, j) 
Endif 

75 Else ! Point is non-dirt, no interpolation! 

!Copy input GLs to output: 
Do ib*l,nb 

im2 (ib, i, j) «iml (ib, i, j) 

Enddo 

cls2(i, j)«cls(i, j) 

- Endif 

20 Enddo !End "i" loop 

Enddo !End "j" loop 

write (6, 500) dcnt,cnt 
500 formatC No. of dirt pels:',i6,'; no. of cleaned pels:',i6) 

If (dent .gt. 0 -and. cnt .eq. 0) then 

write (6,*)' Warning — stable w/ some dirt remaining!!!' 
write (6, *>' [aborting]' 
Call Exit(l) 

Endif 

return 

end 
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13 

c irgb_write_isl: Writes a 3-band image from a Int*2 array- 
Subroutine Irgb_write isl (bname,b, nbx, nb, npixm, nlinm, npix, nlin, 
5 l pixf orm, lunprC) 

Character*50 bname, outf il 

Integer* 2 b (nbx, npixm, nlinm) , npix, nlin, ibuf f (40 96/ 3) 
integer imgout, pstatus, ostatus 
integer nbands, pixf orm, nbx, nb, npixm, nlinm 
integer wstatus, npx, nln, ival 
integer arrtyp, lininc, jtype,pixdira, lunprt 
io integer putdef , opnimg, wrline 

Data pixdim/4096/ 

include ' /local/include/ iopackagef 

arrtyp=IDINT2 
nb an ds=nb 

75 

write (6, 45} npix, nlin, nbands , pixf orm 
45 format <2x, ' w Irgb_write_isl" : : Npix:', 15,'; Nlm:',i5,'; Nb: ,i2, 

1 '; Type:',i2) 

If (lunprt . gt . 0) then 

write (lunprt, 45) npix, nlin, nbands, pixf orm 

c--.Open output image: 

imgout=0 
npx=npix 
nln=nlin 
outfil=bname 

25 pstatus=putdef (imgout, npx, nln, nbands, pixf orm) 

If (pstatus .ne. SYSNRM) then 

write (6,*)' Pstatus error for output — abort' 
stop 

Endif 

ostatus=opnimg (imgout, outf il,l, .false.) !"1" for write-only! 
If (ostatus .ne. SYSNRM) then 
3 0 write (6,*)' Output image not open — abort' 

stop 

Endif 

c Write each line of output image: 

lininc=nlin/5 
jtype=lininc 
35 Do j=l,nlin 

c load 2-d array & write lines: 

Do i=l,npix 

Do ib=l, nbands 

ibuf f (i, ib) =b (ib, i, j) 

Enddo 
Enddo 

40 If (pixform .eq. IDBYTE) then 

Do ib=l, nbands 
Do i=l,npix 

ival^MaxO (ibuff (i,ib) , 0) 
ibuf f (i, ib) *MinO (ival , 255) 

Enddo 

>»c Enddo 

Else if (pixform .eq. IDINT1) then 
Do ib=l, nbands 
Do i*=l,npix 

ival-MaxO (ibuff <i,ib) , -128) 
ibuff (i, ib) «MinO(ival, 127) 

Enddo 

50 Enddo 
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20 


14 

Else if (pixronn .eq. IDINT2) thsn 
continue 

Else if (pixform -eq. IDREA4) then 
continue 

Else 

write (6, 100) pixform 
100 formate Pixf ormo=' , i3, ' not implemented—stop ) 

stop 

Enciif 

line=j -1 _ t 

wstatus=wrline (imgout , line, -1, ibuf f , pixdim, arrtyp) 
If (wstatus .ne. SYSNRM) then 

write (6,*)' Error in wstatus — abort' 

stop 

Endif 

- If (j .ge. jtype) then 
Z write<6,200) j 

c 200 format <10x, ' Line ',i4,' has been written ) 

c jtype=jtype +lininc 

c Endif 

Enddo 'End line loop! 

c. Close image: 

Call els img (imgout) 

return 

end 


c 


25 

Claims 

1. A method for the detection and removal of defects in digital images comprising the steps of: 

a) preprocessing a digital image represented by pixel values to form a feature image in which defect 
30 pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of 
pixel values expected of defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values that are within the range 
of expected defects; 

35 d) forming a residual value for each pixel as a function of the difference between the corresponding 

feature image pixel value and the filtered value of step c); 

e) testing each residual value of step d) to determine if the residual value is within a range of 
. residuals values expected of defects; 

f) forming a map in which each pixel that is within the range of expected defects according to steps 
40 b) and e) is marked as a defect; and 

g) correcting the digital image using the map created in step f). 

2. A method for the detection and removal of defects in digital images comprising the steps of: 

a) preprocessing a digital image represented by pixel values to form a feature image in which defect 
45 pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of 
pixel values expected of defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values that are within the range 
of expected defects; 

so d) forming a residuals image which consists of the difference between the feature image pixel values 

and the filtered values of step c); 

e) testing each pixel value of the residuals image of step d) to determine if the residual value is 
within a range of residuals values expected of defects; 

f) forming a map in which each pixel that is within the range of expected defects according to steps 
55 b) and e) is marked as a defect; and 

g) correcting the digital image using the map created in step f). 

3. A method for the detection and removal of defects in digital images comprising the steps of: 
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a) preprocessing a digital image represented by pixel values to form a feature image in which defect 
pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of 
pixel values expected of defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values; 

d) forming a residual value for each pixel as a function of the difference between the corresponding 
feature image pixel value and the filtered value of step c); 

e) testing each residual value of step d) to determine if the residual value is within a range of 
residuals values expected of defects; 

f) forming a map in which each pixel that is within the range of expected defects according to steps 
b) and e) is marked as a defect; and 

g) correcting the digital image using the map created in step f). 

The method according to Claim 1, 2, or 3, wherein the map formed in step f) is conditioned by 
morphological dilation of the pixels marked as defect. 

The method according to Claim 1, 2, or 3, wherein the map formed in step f) is conditioned by 
morphological erosion of the pixels not marked as defect. 

The method according to Claim 1, 2, or 3, wherein the map of step f) is formed with each pixel marked 
according to the magnitude of its values according to steps b) and e). 

A method for the detection and removal of defects in digital images represented by multiple spatial 
resolutions ordered from lower to higher, comprising the steps of: 

a) preprocessing a lower spatial resolution representation of the digital image by forming a feature 
image in which defect pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of 
pixel values expected of defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values that are within the range 
of expected defects; 

d) forming a residual value for each pixel as a function of the difference between the corresponding 
feature image pixel value and the filtered value of step c); 

e) testing each residual value of step d) to determine if the residual value is within a range of 
residuals values expected of defects; 

f) forming a map in which each pixel that is within the range of expected defects according to steps 
b) and e) is marked as a defect; 

g) correcting the lower resolution representation of the digital image using the map created in step 
0; and 

h) correcting the higher resolution representations of the digital image using the map created in step 
f). 

An apparatus for the detection and removal of defects in digital images comprising: 

means for preprocessing a digital image represented by pixel values to form a feature image in 
which defect pixels have large local contrast from their neighboring non-defect pixels; 

first means for testing the value of each feature image pixel value to determine if the value is within 
a range of pixel values expected of defects; 

filter means for edge-preserving spatial filtering of each of the feature image pixel values; 

means for forming a residual value for each pixel as a function of the difference between the 
corresponding feature image pixel value and the filtered value from said filter means; 

second means for testing each residual value from said forming means to determine if the residual 
value is within a range of residuals values expected of defects; 

means for forming a map in which each pixel that is within the range of expected defects according 
to said first and said second testing means is marked as a defect; and 

means for correcting the digital image as a function of the formed map. 

An apparatus for the detection and removal of defects in digital images represented by multiple spatial 
resolutions ordered from lower to higher, comprising: 

means for preprocessing a lower spatial resolution representation of the digital image by forming a 
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feature image in which defect pixels have large local contrast from their neighboring non-defect pixels; 

first testing means for testing the value of each feature image pixel value to determine if the value 
is within a range of pixel values expected of defects; 

filter means for edge-preserving spatial filtering of each of the feature image pixel values that are 
within the range of expected defects; 

residual means for forming a residual value for each pixel as a function of the difference between 
the corresponding feature image pixel value and the filtered value of said filter means; 

second testing means for testing each residual value from said residual means to determine if the 
residual value is within a range of residuals values expected of defects; 

mapping means for forming a map in which each pixel that is within the range of expected defects 
from said first and said second testing means is marked as a defect; 

means for correcting the lower resolution representation of the digital image using the map created 
by said mapping means; and 

means for correcting the higher resolution representations of the digital image using the map 
created by said mapping means. 

10. A method for the detection and removal of local defects in digital images comprising the steps of: 

a) generating a feature image from an original digital image; 

b) creating an EPS-residuals image using the feature image; 

c) testing and mapping image pixels based upon the values of the EPS residuals image and the 
pixel values of the feature image; 

d) modifying the mapped image as a function of mapped values of neighboring pixels; 

e) cleaning of defect pixels by an EPS-substitution process or by a greyscale erosion process; and 

f) cleaning of defect pixels in different spatial resolution representations of the original image if 
different spatial resolution representations are present. 

11. A method for the detection and removal of defects in digital images comprising the steps of: 

a) preprocessing a digital image to form a feature image; and 

b) testing each feature image to determine if a value is within a range of pixel values expected of 
defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values that are within the range 
of expected defects. 
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